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A  . 


ABSTRACT 

A  two-dimensional,  two-passage  simulation  of  the  relative  flow  through  a  transonic 
fan  at  M  =  1.4  was  designed  with  a  view  to  providing  an  apparatus  in  which  to  assess  the 
effectiveness  of  passive  vortex  generator  techniques  in  alleviating  shock-boundary  layer 
interaction  effects.  The  design  of  the  model  and  the  results  of  six  tests  in  the  transonic 
cascade  blowdown  wind  tunnel  are  described.  Schlieren  photographs  of  the  shock 
structure  were  obtained  at  back  pressures  lower  than  the  design  value.  The  back-pressure 
control  valve  was  identified  as  being  critical  to  completing  the  experimental  simulation. 
The  flow  through  the  cascade  geometry  was  computed  at  design  pressure  ratio  using  an 
Euler  code.  Modifications  to  the  grid  are  recommended  before  thin-layer  Navier-Stokes 
calculations  are  performed. 
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I.    INTRODUCTION 

A.      SHOCK-BOUNDARY  LAYER  INTERACTION 

The  motivation  for  the  present  study  is  as  introduced  in  Ref.  1 : 

Outboard  of  approximately  50-60%  span  on  current  fan  designs,  the  peak  relative 
Mach  Number  reaches  1.3.  There  is  a  normal  shock  wave  that  stands  in  front  of  the 
passage  that  is  of  sufficient  strength  to  separate  the  boundary  layer  on  the  suction  side, 
as  shown  in  Fig.  1.  Depending  on  the  strength  of  the  shock-boundary  layer  interaction, 
the  boundary  layer  may  or  may  not  reattach  itself  to  the  suction  side.  As  a  result  of  this 
separation,  there  is  little  or  no  pressure  recovery  downstream  of  the  shock  and  there  are 
high  total  pressure  losses.  If  this  separation  is  able  to  be  suppressed,  which  would 
improve  the  boundary  layer  characteristics  downstream  of  the  shock,  the  stage  pressure 
rise  would  increase  and  the  fan  efficiency  would  improve.  In  addition,  the  fan  rotation 
speed  could  increase,  allowing  a  reduction  in  size  and  weight  and  more  efficient  turbine 
speeds. 

A  diagram  of  the  basic  flow  features  of  the  shock-boundary  layer  interaction  are 
shown  at  the  bottom  of  Fig.  1.  The  step  increase  across  the  shock  cannot  be  handled  by 
the  subsonic  portion  of  the  boundary  layer,  therefore,  this  pressure  rise  is  partially 
transmitted  upstream  through  the  subsonic  region.  This  results  in  the  divergence  of  the 
stream  lines  and  formation  of  compression  waves  upstream  of  the  shock.    If  the  shock 


pressure  rise  is  sufficiently  strong,  the  boundary  layer  will  separate  and  cause  a  more 
severe  thickening  of  the  upstream  boundary  layer.  This  will  cause  the  coalescence  of  the 
compression  waves  into  the  leading  shock.  In  order  for  the  flow  to  turn  back  towards  the 
axial  direction,  a  rear  shock  is  formed.  The  leading  and  rear  shock  combination  is  called 
the  lambda  foot.  The  total  pressure  loss  through  the  lambda  foot  is  lower  than  the  loss 
across  the  normal  shock,  resulting  in  a  free  shear  layer  originating  at  the  shock  triple 
point.  The  undisturbed  boundary  layer  thickness  at  the  shock  location,  80,  is  the 
normalizing  length  used  in  describing  the  geometry  of  the  interaction. 

Recent  experiments  have  examined  several  concepts  which  have  been  proposed  for 
suppressing  the  induced  separation  caused  by  the  shock-boundary  layer  interaction.  Three 
methods  investigated  at  the  United  Technologies  Research  Center  (UTRC)  include  low- 
profile  vortex  generators,  a  porous  wall  with  a  passive  cavity  and  surface  contouring  [Ref. 
1].  Additional  methods  examined  at  the  NASA  Langley  Research  Center  include  large 
eddy  breakup  Devices  (LEBU)  at  small  angle  of  attack,  spanwise  cylinders  and  vortex 
generator  jets  [Ref.  2].  The  most  promising  and  adaptable  for  modern  transonic  cascades 
are  the  vortex  generators  and  the  vortex  generator  jets. 

B.      VORTEX  GENERATING  DEVICES 

One  proposed  method  for  supressing  or  eliminating  the  shock-induced  separation 
is  by  the  replenishment  of  the  boundary  layer  momentum  near  the  surface  [Ref.  2]. 
Replenishment  can  be  effected  by  a  redirection  of  the  momentum  in  the  free  stream  or 
outer  region  of  the  boundary  layer  into  the  near-wall  flow.   A  mechanism  for  effecting 


the  momentum  transfer  is  via  the  introduction  of  streamwise  vorticity.  Such  a  momentum 
exchange  in  the  boundary  layer  can  enable  the  layer  to  adjust  to  the  sudden  pressure  rise 
without  separation  [Ref.  1].  Various  techniques  can  be  used  to  introduce  streamwise 
vorticity.  Two  techniques,  vortex  generators  and  vortex  generator  jets,  are  described  in 
the  following  sections. 

1.      Vortex  Generators 

Passive  vortex  generators  with  a  height  on  the  order  of  the  boundary  layer  thickness 
have  long  been  known  to  increase  the  mixing  between  external  streams  and  the  boundary 
layer.  Such  vortex  generators,  producing  relatively  large  streamwise  trailing  vortices, 
have  been  used  effectively  to  delay  boundary  separation,  to  enhance  aircraft  wing  lift  and 
to  avoid  or  delay  separation  in  subsonic  diffusers  [Ref.  2].  However,  because  of  the 
fullness  of  the  mean  velocity  profile  of  high  Reynolds  Number  turbulent  boundary  layers, 
a  reduction  of  the  height  of  the  vortex  generators  to  a  fraction  (25-50%)  of  the  boundary 
layer  height  has  been  proposed  [Ref.  2].  The  small  profile  offers  smaller  perturbations 
to  the  supersonic  freestream.  Also,  they  are  rugged  in  design,  capable  of  withstanding 
centrifugal  loads  and  other  harsh  conditions.  Lastly,  they  have  demonstrated  a  working 
flow  incidence  of  ±5°  [Ref.  1].  The  low-profile  vortex  generators,  shown  in  Fig.  2,  were 
tested  in  a  flow  in  which  the  boundary  layer  was  approximately  1.28  in.  thick.  They 
include  the  doublet  and  the  wishbone  [Ref.  2].  Other  designs  include  the  asymmetric 
candleflame  and  the  singlet  (a  doublet  with  the  second  wedge  removed)  [Ref.  1]. 

UTRC  investigated  the  effects  of  low-profile  vortex  generators  on  suppressing 
shock-induced  separation  in  a  slowly  diverging  circular  duct.   Details  of  the  experiment 


are  found  in  Ref.  1.  The  shock  position  within  the  duct  was  controlled  using  a  translating 
plug  valve  following  the  test  section.  Wall  pressure  distributions  through  the  interaction 
region  were  obtained  and  surface  flow  structure  was  observed  using  an  oil  flow  technique. 

The  results  showed  that  the  pressure  rise  increased  towards  the  ideal  value  when 
vortex  generators  were  attached  upstream,  and  the  region  of  the  pressure  rise  was 
shortened.  Thus  the  shock-induced  separation  was  substantially  reduced.  Depending  on 
the  configuration  of  vortex  generators,  the  shock-boundary  layer  interaction  length  could 
be  shortened  by  approximately  50%.  Also,  the  displacement  and  momentum  boundary 
layer  thicknesses  downstream  were  substantially  reduced.  A  sketch  which  illustrates 
changes  in  the  pressure  recovery  obtained  with  vortex  generators  is  shown  in  Fig.  3  in 
comparison  with  the  ideal  pressure  recovery  and  the  baseline  pressure  recovery. 

Of  the  remaining  methods  of  suppressing  shock-induced  separation  investigated  in 
Ref.  1,  the  porous  wall  with  the  passive  cavity  did  reduce  the  shock  loss,  but  the 
boundary  layer  development  was  aggravated.  Also,  significant  mixing  losses  occured 
downstream.  The  use  of  surface  contouring  moderately  reduced  the  boundary  layer 
thickness,  but  did  not  affect  the  size  of  the  separation  region. 

2.      Vortex  Generator  Jets 

The  use  of  vortex  generator  jets  was  first  proposed  by  Wallis  [Ref.  3]  in  Australia. 
The  vortex  jets  are  usually,  but  not  as  intended  here,  an  active  method  where  jets  of  air 
are  blown  through  holes  in  the  solid  surface.  The  holes  are  inclined  at  an  angle  to  the 
surface  and  skewed  with  respect  to  the  fireestream.  A  schematic  of  the  arrangement 
(taken  from  Ref.  2)  is  shown  in  Fig.  4.    The  arrangement  of  the  jets  on  the  surface  is 


much  like  the  arrangement  of  solid  vortex  generators.  The  relatively  weak  counter- 
rotating  vortices  that  are  formed  within  the  jet  fluid  in  crossflow  are  replaced  by  a 
stronger  streamwise  vortex  trail  in  the  ambient  fluid  downstream  of  injection  and  close 
to  the  surface  [Ref.  2]. 

Vortex  generator  jets  have  proven  effective  in  reducing  the  extent  of  the  turbulent 
boundary  layer  separation  resulting  from  adverse  pressure  gradients.  However,  the  status 
is  summarized  in  the  following  points  stated  in  Ref.  4: 

•  The  feasibility  of  the  technology  is  not  firmly  established 

•  It  is  possible  to  adapt  the  arrangement  passively  to  fans  because  of  the  cascade 
arrangment 

•  There  is  a  negligible  drag  penalty  of  the  vortex  generator  jets  when  compared  to  the 
vortex  generators  when  the  jets  are  off. 

The  effect  of  skewed  jets  and  their  relative  angles  (6)  to  the  freestream  on 
separation  in  a  subsonic  diffusing  flow  were  investigated  in  Ref.  4.  The  jets  were  skewed 
45°  to  the  vertical  and  6  was  varied  from  0°  (downstream)  to  180°  (upstream).  The 
following  conclusions  were  reached. 


•  Substantial  reduction  of  the  size  of  the  stall  region  size  occurred  for  jet  velocity 
ratios  >  0.8. 

•  Jets  pointing  directly  upstream  were  ineffective  due  to  weak  vortex  formation. 

•  Jet  arrays  that  cause  counter-rotating  vortex  pairs  can  cause  significant  spanwise 
variations. 

•  Spanwise  average  reattachment  location  seems  to  be  more  strongly  affected  than  the 
detachment  location 


C.      2-D  FAN  PASSAGE  SIMULATION 

The  promising  results  of  the  UTRC  experiment  suggested  the  need  for  an  evaluation 
of  the  concept  in  an  actual  transonic  fan  configuration.  A  full-scale  simulation  would  be 
desirable,  but  such  an  experiment  has  the  potential  problem  of  not  being  definitive 
because  of  extreme  difficulty  in  measurement  and,  therefore,  in  adequately  assessing  the 
results.  Hence,  a  steady  2-D  stream- surface  simulation  was  undertaken  in  the  present 
study.  The  combination  of  an  experimental  and  computational  simulation  was  the 
approach  taken  in  an  attempt  to  establish  a  representative  baseline  flow  geometry  in 
which  to  verify  the  effectiveness  of  the  proposed  vortex  generators  in  alleviating  the 
shock-induced  separation. 

The  geometry  for  the  2-D  experiment  was  a  simulation  of  the  relative  flow  on  an 
advanced  fan  rotor  at  approximately  63%  of  the  span.  The  2-D  model  used  in  the 
experiment  was  based  on  the  stream  surface  geometry,  but  the  blade  profile  was 
approximated  (closely)  as  a  wedge  arc  for  ease  in  manufacture.  This  was  a  logical  step 
since  steamline  contraction  could  not  be  simulated  in  the  experiment.  The  geometry  of 
the  2-D  simulation  is  shown  in  Fig.  5. 

In  reporting  the  generation  of  the  baseline  flow  geometry  in  the  present  document, 
Section  II  describes  the  experimental  simulation  and  Section  HI  describes  the  numerical 
simulation.  As  originally  conceived,  from  the  experimental  simulation,  measurements  of 
surface  pressures  and  schlieren  photographs  of  the  shock-boundary  layer  interaction  were 
to  be  obtained.    Measured  boundary  pressures  were  then  to  be  used  as  input  boundary 


conditions  for  the  numerical  simulation,  and  a  comparison  could  be  made  of  both  shock- 
boundary  layer  interaction  structure  and  distribution  of  surface  pressure.  However,  as 
described  in  the  results  given  in  Section  IV,  while  the  experimental  apparatus  was 
relatively  successful,  it  could  not  be  operated  at  sufficiently  high  back-pressures  (>2x 
upstream)  to  produce  a  proper  simulation  of  the  fan  shock  structure  at  design  conditions. 
Also  given  in  Section  IV  are  the  results  of  a  computation  of  the  flow  for  the  design 
pressure  ratio  of  2.28  using  an  Euler  code.  It  was  noted  that,  in  the  absence  of  viscosity 
in  the  code  calculations  at  design  pressure  ratio,  and  with  too  low  a  back-pressure  in  the 
experiment,  both  computation  and  schlieren  results  showed  qualitative  agreement  for 
shock  structure  at  the  entrance  to  the  passage.  It  was  concluded,  in  Section  V,  that 
improvement  (through  a  refinement  in  the  design  of  the  valve)  in  the  control  of  the  back- 
pressure is  required  to  properly  simulate  fan  conditions.  Grid  refinement  is  required  to 
continue  the  computational  simulation. 

Details  of  the  experiment  are  given  in  Appendices  A  -  C  and  details  of  the 
numerical  approach  and  the  grid  generation  code  are  given  in  Appendices  D  -  F. 


II.   EXPERIMENTAL  SIMULATION 

A.      TRANSONIC  CASCADE  WIND  TUNNEL 

1.       Wind  Tunnel  Description 

The  original  specifications,  design  and  construction  of  the  Transonic  Cascade  Wind 
Tunnel  are  described  in  Ref.  5.  Subsequent  modifications  to  the  tunnel  are  described  in 
Ref.  6. 

The  apparatus  is  a  blowdown  device  with  a  design  Mach  number  of  1.4  and  a 
design  total  pressure  of  50  psia.  A  schematic  of  the  Transonic  Cascade  Wind  Tunnel  is 
shown  in  Fig.  6.  The  tunnel  was  relocated  (to  Bldg.  216  at  NPS)  and  operated  with  a 
new  300  psia  air  supply  system  in  the  course  of  the  present  work.  A  photograph  of  the 
right  hand  side  of  the  wind  tunnel  is  shown  in  Fig.  7,  and  a  photograph  of  the  left  hand 
side  is  shown  in  Fig.  8.  A  tubular  flow  straightener  is  used  to  remove  the  swirl  induced 
by  the  upstream  pipe  bends  and  the  control  valve.  The  settling  chamber  is  composed  of 
a  flat,  perforated  plate  diffuser  followed  by  a  screen  to  reduce  the  turbulence  level  and 
promote  flow  uniformity.  There  is  a  round  (10  in.  diameter)  to  rectangular  (2.0  in  x  4.5 
in)  transition  section  from  the  settling  chamber  to  the  supersonic  nozzle.  The  tunnel  uses 
a  2-dimensional  converging-diverging  nozzle  to  accelerate  to  M  =  1.4.  A  photograph  of 
the  interior  of  the  test  section  (without  blades)  is  shown  in  Fig.  9.  Aft  of  the  test  section 
is  a  back  pressure  valve  that  is  used  to  adjust  the  shock  location  within  the  test  model. 
The  back  pressure  valve  is  a  "ramp  and  drum"  throttle  valve  with  the  ramp  actuated  after 
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the  flow  "starts"  by  a  pneumatic  cylinder.  A  schematic  of  the  back  pressure  valve  is 
shown  in  Fig.  10,  and  a  photograph  of  the  exterior  of  the  valve  is  shown  in  Fig.  11.  The 
exhaust  is  ducted  to  a  sound  suppressor  outside  of  the  building.  The  tunnel  start-up  is 
controlled  by  a  pneumatic  control  valve.  A  photograph  of  the  control  valve  is  shown  in 
Fig.  12,  and  a  photograph  of  the  control  panel  is  shown  in  Fig.  13. 

2.       Data  Acquisition  System 

The  data  acquisition  system  used  in  the  experiment  is  described  in  detail  in  Ref.  7. 
A  schematic  of  the  system  is  shown  in  Fig.  14,  and  a  view  of  the  system  is  shown  in  Fig. 
15.  The  BASIC  code  used  to  record  the  tunnel  reference  pressure  measurements  (using 
an  uninstrumented  bottom  blade)  is  given  in  Appendix  B,  Fig.  Bl.  The  Scanivalve  port 
number  and  corresponding  pressure  location  were  as  follows: 

•  Port  1  -  Atmospheric  pressure 

•  Port  2  -  Calibration  pressure  (50  psia) 

•  Port  4  -  Plenum  pressure 

•  Port  6  -  Upstream  static  pressure 

•  Port  8  -  Downstream  static  pressure 

Intervening  ports  were  vented  to  the  atmosphere  to  verify  that  the  measured  pressures  had 
time  to  become  steady  state. 
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3.       Optical  System 

A  schematic  of  the  optical  system  is  shown  in  Fig.  16.  The  system  included  a 
continuous  light  source  schlieren  system  and  a  spark  source  schlieren  system.  The 
continuous  light  source  was  a  100  W  mercury  vapor  arc  lamp,  and  the  spark  source  was 
a  200  -  300  nsec  spark  discharge.  The  parabolic  lens  provided  a  parallel  light  beam 
through  the  test  section.  The  parabolic  mirror,  mirror  2  in  Fig.  16,  on  the  far  side  of  the 
test  section  focused  the  parallel  beam  onto  the  knife  edge  located  in  front  of  the  camera 
box  unit.  At  the  rear  of  the  camera  box  was  a  viewing  screen  behind  a  shutter  capable 
of  automatic  or  manual  operation.  The  minimum  exposure  time  was  one-thousandth  of 
a  second.  A  Polaroid  film  holder  could  be  inserted  into  the  camera  box  when  needed. 
In  the  present  experiment  black-and-white  Polaroid  Type  52  film  was  used.  A  front 
surface  flat  mirror,  mirror  1  in  Fig.  10,  was  inserted  between  the  continuous  light  source 
and  the  parabolic  lens  when  the  spark  source  was  to  be  operated.  Mirror  one  was 
positioned  and  aligned  using  a  highly  focused  flashlight  Both  the  continuous  light  source 
and  the  spark  source  were  filtered  to  provide  high  resolution  schlieren  photographs.  A 
view  of  the  left  hand  side  of  the  optical  system  is  shown  in  Fig.  17. 

B.      TEST  SECTION  DESIGN 

1.       Design  of  the  Blading 

Based  on  stream-surface  conditions  through  the  3-D  fan  blade,  the  2-D  experiment 
was  designed  to  retain  the  approximate  geometry  while  simplifying  curvatures  for  ease 
of  manufacture.     It  was  found  that  the  stream-surface  fan  blade  profile  could  be 
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approximated  by  a  wedge-arc  blade  shape  having  a  flat  pressure  side  and  equal  leading 
and  trailing  edge  radii  as  listed  in  Fig.  5.  The  leading  and  trailing  edge  radii  were  0.015 
in,  and  the  chord  of  the  blade  was  6.000  in.  The  suction  surface  starting  at  the  leading 
edge  was  inclined  to  the  pressure  surface  at  3.5°  to  a  point  2.85  in.  along  the  pressure 
surface.  The  suction  surface  wedge  was  blended  to  the  trailing  edge  radius  by  an  arc  of 
radius  13.53  in. 

The  test  section  geometry  (oriented  to  allow  the  flow  direction  to  be  horizontal)  was 
composed  of  three  partial  blades  defining  two  blade  passages,  as  shown  in  Fig  18.  The 
top  passage  was  to  be  the  control  passage  while  the  bottom  passage  was  to  be  used  to 
study  the  effects  of  the  vortex  generators.  For  the  experiment  with  the  vortex  generator 
jets,  the  bottom  blade  was  for  control,  while  holes  were  to  be  drilled  in  the  center  blade. 
A  photograph  of  the  model  blade  is  shown  in  Fig.  19,  while  the  machine  drawing  of  the 
model  blade  is  given  in  Appendix  A,  Fig.  A 10.  Figures  20  and  21  show  views  of  the 
upper  and  lower  partial  blades  respectively,  while  the  machine  drawings  of  the  blades  are 
given  in  Appendix  A,  Figs.  A9  and  All  respectively.  Appendix  A,  Fig.  A 12,  details  the 
location  of  the  pressure  taps  for  the  instrumented  lower  blade.  Figure  22  shows  all  three 
blades  attached  to  the  inner  plate. 

2.       Design  of  the  Side  Pieces,  Inner  Plates  and  Windows 

The  design  of  the  test  section  was  driven  by  three  major  considerations:  nozzle  wall 
boundary  layer  diversion  from  the  nozzle,  wave  cancellation  from  the  leading  edges  of 
the  blades  and  optical  accessibility.  Boundary  layer  diversion  for  the  top  and  bottom  of 
the  test  section  was  accomplished  through  bleed  channels  above  the  top  blade  and  below 
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the  bottom  blade,  using  existing  ducting  [Ref.  5].  Side  wall  boundary  layer  diversion  was 
provided  by  designing  an  inner  plate  which  began  near  the  leading  edges  of  the  blades. 
Figure  19  shows  a  photograph  of  the  inner  plate  and  center  blade,  while  Appendix  A,  Fig. 
A2,  shows  the  machine  drawing  of  the  inner  plate.  The  exhaust  for  the  diverted  side  wall 
boundary  layer  was  through  a  slot  in  the  side  piece.  Figure  23  shows  a  photograph  of  the 
side  piece  and  Appendix  A,  Fig.  Al,  shows  the  machine  drawing  of  the  side  piece. 
Appendix  A,  Fig.  A4,  details  the  design  of  the  side  piece  exhaust  slot.  In  addition,  the 
side  piece  and  inner  plate  incorporated  a  window  for  optical  measurements.  The  window 
was  made  of  3/4  in.  Plexiglas,  and  was  designed  to  allow  viewing  access  to  the 
anticipated  shock-boundary  interaction  regions  in  both  passages.  Preliminary  tests 
indicated  that  the  thick  Plexiglas  would  give  acceptable  schlieren  quality.  Lastly,  the 
cancellation  of  weak  waves  from  the  leading  edge  was  effected  by  adjusting  bleed  through 
a  porous  wall  on  the  top  of  the  test  section  forward  of  the  top  blade  [Ref.  6  and  8]. 

C.      EXPERIMENTAL  RESULTS 

Six  test  runs  were  made.  The  first  experiment,  Run  001,  was  with  the  back  pressure 
valve  initially  completely  open.  The  flow  angle  of  incidence  (to  the  suction  surface)  was 
set  at  1.15°.  Figure  24  shows  the  continuous  light  schlieren  photograph  which  was 
obtained.  The  photograph  shows  an  oblique  shock  leading  from  the  leading  edge  of  the 
middle  blade  to  the  suction  surface  of  the  bottom  blade.  The  multiple  waves  on  the  top 
blade  are  from  the  porous  wall,  which  was  not  adjusted  to  reduce  reflections.  The 
parallel  weaker  oblique  shock  in  the  lower  passage  is  believed  to  originate  at  the 
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intersection  of  the  leading  edge  of  the  supporting  tab  on  the  middle  blade  and  the  inner 
frame.  Once  the  photograph  was  taken  the  back  pressure  valve  was  slowly  closed  in  an 
effort  to  bring  the  oblique  shock  forward.  Figure  25  shows  the  continuous  light  schlieren 
photograph  when  the  back  pressure  valve  was  rotated  to  fully  close  the  gap  between  the 
eccentrically  mounted  cylinder  and  the  (fully  raised)  ramp.  The  oblique  shock  is  seen  to 
have  moved  forward  towards  the  leading  edge  in  both  passages.  The  shock  appears  to 
be  near  normal  impinging  on  the  boundary  layer.  Tables  IA  and  IB  lists  the  pressure 
recordings  for  Run  001. 

The  experimental  pressure  ratio  p2/p!  was  found  to  be  0.765  with  the  back  pressure 
valve  open.  With  the  back  pressure  valve  fully  closed,  p2/pi  =  1.39.  Calculations  for 
the  fan  geometry  indicated  that  the  pressure  ratio  should  be  2.28  on  design  when  a  strong 
normal  shock  would  exist  across  the  blade  passage.  On  inspection,  the  throat  of  the  valve 
in  the  fully  closed  position  when  added  to  the  spaces  around  the  sides  of  the  ramp  gave 
an  area  which  was  too  large  to  achieve  the  desired  throttling. 

Run  002  was  made  after  the  actuator  was  adjusted  to  allow  the  ramp  and  drum  to 
fully  close  when  the  drum  was  actuated.  However,  the  pressure  ratio  p2/p,  obtained  in 
Run  002  was  1.26.  No  photographs  were  taken  of  this  run.  Run  003  was  made  with  the 
actuation  pressure  increased  to  600  psia.  Again,  inadequate  pressure  ratios  were  obtained. 
The  schlieren  photographs  and  pressure  measurements  (not  shown)  were  very  similar  to 
those  recorded  in  Run  001. 

Run  004,  with  the  ramp  actuator  operated  at  750  psia,  showed  no  improvement. 
Run  005  was  run  with  the  back-pressure  valve  mechanically  fixed  in  a  "closed"  position 
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using  a  pipe  sleeve  around  the  actuator  rod.  Of  concern  was  the  possibility  of  the  flow 
not  starting.  Although  the  flow  was  started  successfully  through  the  model,  it  was  found 
that  the  ramp  was  1/16  -  1/8  in.  lower  than  was  required  to  fully  close,  and,  again,  the 
required  control  was  inadequate.  For  Run  006,  the  ramp  was  wedged  in  the  fully  up 
position  and  the  spaces  between  the  ramp  and  duct  walls  were  sealed  using  composite 
material  plates.  During  this  test,  a  maximum  back-pressure  ratio  p2/p,  =  1.98  was 
obtained  before  the  wedge  holding  the  pipe  in  place  and  one  composite  plate  gave  way. 
As  the  back-pressure  approached  1.98,  the  shock  location  moved  toward  the  leading  edge 
of  the  blade;  however,  the  motion  was  very  unsteady  and  very  sensitive  to  the  back- 
pressure valve  position.  Table  IC  lists  the  pressures  recorded  in  Run  006.  No  schlieren 
photograph  was  attempted  in  Run  006  in  view  of  the  need  to  observe  the  shock  movement 
while  experimenting  with  the  valve  position. 
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in.   NUMERICAL  SIMULATION 

A.  COMPUTATIONAL  SCHEME 

The  computational  code  used  in  the  present  work  was  based  on  the  thin-layer 
Navier-Stokes  equations.  The  scheme  used  a  body-fitted  coordinate  system,  and  was  a 
full  upwind  algorithm  based  on  Roe  flux-difference  splitting.  The  Baldwin-Lomax 
turbulence  model  was  used  except  near  shock  impingement.  The  Johnson-King 
turbulence  model  was  used  to  model  the  shock-induced  separation.  Lastly,  the  Osher- 
Chakravarthy  TVD  scheme  for  flux-limiting  was  used.  Details  of  the  development  of  the 
numerical  code  and  background  theory  are  found  in  Appendix  D. 

B.  GRID  GENERATION 

The  grid  for  the  numerical  scheme  was  generated  using  the  GRAPE  grid  generation 
code,  described  in  Appendix  E.  Two  different  grids  were  generated:  a  channel  passage 
C-grid,  to  be  used  to  carry  out  a  preliminary  calculation  of  the  flow  over  a  single  blade 
between  parallel  walls,  and  a  cascade  C-grid,  to  be  used  to  compute  the  flow  through  the 
simulated  fan  passage. 

The  channel  passage  C-grid  was  an  199  x  64  grid  suitable  for  Euler  calculations. 
The  input  file  for  use  with  GRAPE  is  shown  in  Appendix  F,  Fig.  Fl.  Figure  26  shows 
the  channel  passage  C-grid. 
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The  cascade  C-grid  was  generated  using  a  modification  made  to  GRAPE  which  is 
outlined  in  Appendix  E.  The  input  file  for  GRAPE  is  shown  in  Appendix  F,  Fig.  F2. 
The  grid  size  was  169  x  40.  The  cascade  feature  was  used  to  obtain  the  necessary  flow 
periodicity  in  the  developed  grid.  Figure  27  shows  the  complete  cascade  C-grid,  while 
Fig.  28  shows  an  enlarged  view  of  the  leading  edge  region.  Figures  29  and  30  show  the 
details  near  the  leading  and  trailing  edges  respectively.  Figure  31  shows  two  cascade  C- 
grids  plotted  simultaneously  with  one  grid  at  an  offset  equal  to  the  blade  spacing.  This 
figure  demonstrates  the  periodicity  of  the  grid  and  the  accuracy  of  the  geometric 
specifications  such  as  the  solidity.  Figure  32  shows  a  detailed  view  of  the  boundary 
between  the  two  grids. 

C.      COMPUTED  RESULTS 

The  Mach  number  contours  for  the  Euler  solution  obtained  for  the  flow  over  the 
single  blade  in  a  channel  using  the  channel  passage  C-grid  are  shown  in  Fig.  33.  The 
solution  shown  in  Fig.  33  is  not  a  steady  state  condition  but  rather  a  "snapshot"  of  the 
flow  during  a  transient  in  which  the  pressure  side  shock  was  propagating  towards  the 
leading  edge. 

The  Mach  number  contours  for  the  Euler  solution  for  the  flow  through  the  cascade 
C-grid  are  shown  in  Fig.  34.  The  angle  of  incidence  of  the  flow  to  the  suction  surface 
was  0°  and  pj/pj  =  2.28.  The  contours  show  a  pressure  side  leading  edge  oblique  shock 
and  a  reflection  from  the  suction  surface,  followed  by  a  strong,  near-normal  passage 
shock.    Figure  35  shows  the  pressure  contours.    On  both  sets  of  contours,  the  leading 
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edge  disturbance  does  not  continue  to  propagate  smoothly  across  grid  boundaries  because 
of  a  discontinuity  encountered  there  in  the  grid  spacing  (see  Fig.  32).  The  calculated 
blade  surface  Cp  distribution  is  shown  in  Fig.  36.  This  figure  also  demonstrates  the 
effectiveness  of  flux-limiters  in  suppressing  the  oscillations  that  develop  around  the 
discontinuities. 
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IV.   DISCUSSION  OF  RESULTS 

A.      EXPERIMENTAL 

The  oblique  shock  patterns  within  the  passage,  shown  in  Figs.  24  and  25,  were 
observed  to  be  qualitatively  similar  to  the  computed  shock  pattern  shown  in  Fig.  34.  It 
is  noted,  however,  that  the  flow  incidence  to  the  suction  surface  in  the  experiment  was 
1.15°  whereas  it  was  0°  in  the  computation.  The  shock  should  therefore  have  been 
stronger  in  the  experiment  than  in  the  computation,  and  this  certainly  was  the  case. 

In  the  present  work,  the  difficulties  with  the  back-pressure  valve  precluded  obtaining 
a  pressure  ratio  p2/p!  =  2.28.  However,  considerable  experience  was  gained  in  the 
attempts  which  were  made  to  achieve  this  ratio.  First,  it  was  found  that  the  flow  could 
be  started  through  the  test  section  with  the  ramp  fully  up.  This  means  that  the  ramp  can 
be  fixed  rigidly  to  the  duct  wall  and  need  not  be  supported  by  a  pneumatic  actuator.  The 
rigid  positioning  will  allow  the  sides  of  the  ramp  to  be  sealed  completely  and  will 
eliminate  a  potential  source  of  unsteadiness  in  the  valve  throat  area. 

For  the  brief  period  during  Run  006  that  the  back-pressure  was  held  at  a  value 
approaching  twice  the  inlet  pressure,  the  flow  was  highly  unsteady.  Unfortunately,  the 
inlet  pressure  was  dropping  as  the  supply  pressure  fell,  and  it  was  not  possible  to 
experiment  further.  It  was  noted  on  the  schlieren  screen  that  the  lower  passage  unstarted 
while  the  upper  passage  remained  started,  raising  the  possibility  that  the  achievement  of 
similar  flows  in  the  two  passages  may  not  be  automatic. 
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B.      NUMERICAL 

The  Euler  solution  using  the  cascade  C-grid  gave  useful  qualitative  information  on 
the  structure  of  the  flow  field,  which  was  supported  by  the  experimental  data.  It  is 
expected  that  at  the  same  pressure  ratio,  the  viscous  solution  should  move  the  shock  in 
the  passage  due  to  losses  and  the  displacement  effects  of  the  boundary  layer.  The  flow 
structure  might  still  not  correspond  to  what  is  expected  to  occur  in  the  fan,  wherein  the 
normal  shock  sits  at  the  leading  edge  in  an  "unstarted"  structure.  It  is  also  possible  that 
the  computed  structure  is  the  result  of  the  initial  conditions  imposed  on  the  unsteady 
solution.  What  is  found  is  a  solution  that  requires  a  "starting"  process  that  is  not  present 
physically  in  the  fan  application.  What  would  be  required  to  obtain  the  "unstarted" 
solution  would  be  to  first  increase  the  pressure  ratio  to  expel  the  shock,  and  then  decrease 
the  pressure  ratio  to  the  design  value.  It  is  interesting  to  note  that  this  is  exactly  the 
procedure  which  must  be  followed  in  the  experiment,  which  is  why  the  back  pressure 
control  valve  is  such  a  critical  component  in  the  test  apparatus. 

Grid  refinement  is  also  necessary  for  the  intended  application  with  a  viscous  flow 
solver.  The  boundary  layer  will  need  to  be  resolved  adequately  with  an  appropriate 
number  of  grid  points.  Also,  the  orthogonality  between  the  upper  and  lower  boundaries 
of  the  cascade  will  have  to  be  enforced  more  rigidly  (Fig.  32).  Lastly,  the  grid  spacing 
will  need  to  be  more  uniform  to  ensure  continuity  across  boundaries  (Fig.  32).  Continuity 
in  grid  spacing  is  expected  to  yield  more  accurate  propagation  of  shock  discontinuities. 
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It  is  noted,  however,  that  the  validity  of  the  turbulence  model  chosen  has  been 
demonstrated  (Appendix  D,  Fig.  Dl),  and  a  solution  for  an  airfoil  has  been  obtained  using 
the  viscous  solver  (Appendix  D,  Fig.  D2). 
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V.   CONCLUSIONS  AND  RECOMMENDATIONS 

In  the  present  study,  both  experimental  and  numerical  simulations  of  the  flow  on 
a  stream  surface  through  a  transonic  fan  passage  were  attempted  with  a  view  to  providing 
a  baseline  geometry  in  which  the  effectiveness  of  upstream  vortex  generating  devices  in 
alleviating  effects  of  shock-boundary  layer  interaction  could  be  examined.  From  the 
experiment,  the  following  were  concluded: 


•  The  design  of  the  test  section  geometry  was  successful.  No  structural  problems 
were  encountered  in  six  test  runs.  Boundary  layer  diversion  slots  appeared  to  start 
by  natural  aspiration  to  the  atmosphere.  Schlieren  photographs  were  obtained 
successfully  through  thick  Plexiglas  windows. 

•  The  upstream  pressure  control  valve  was  fully  satisfactory,  and  the  transonic 
cascade  wind  tunnel  was  made  operative  in  its  new  location. 

•  The  downstream  pressure  control  valve  was  found  not  to  be  satisfactory  as 
originally  installed  and  was  shown  to  be  critical  to  the  eventual  success  of  the 
experiment.  It  was  shown  that  the  test  section  flow  would  establish  with  the  valve 
ramp  fully  up,  so  that  the  ramp  could  be  made  both  rigid  and  sealed.  The  eccentric 
aluminum  cylinder  component  did  not  provide  adequate  control  because  of  internal 
leakage. 

•  No  problem  other  than  the  back-pressure  control  valve  was  found  which  would 
prevent  a  successful  experiment. 


From  the  initial  steps  taken  toward  a  numerical  simulation,  the  following  were 
concluded: 
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•  A  suitable  grid  for  computing  the  cascade  flow  using  an  Euler  solver  was  generated 
using  GRAPE  code  as  modified  by  R.  Chima  (Appendix  E).  A  suitable  grid  for 
viscous  calculations  near  this  blade  surface  can  be  obtained  using  this  code. 
Discontinuities  in  grid  spacing  across  periodic  boundaries  given  by  the  code  can 
possibly  be  removed  using  the  DSOBI  parameter  (Appendix  E). 

•  The  solution  obtained  for  the  experimental  cascade  at  the  design  pressure  ratio 
using  an  Euler  solver  showed  a  normal  shock  wave  well  inside  the  passage.  Since 
the  inclusion  of  viscous  displacement  effects  should  move  the  shock  downstream, 
this  solution  may  be  the  result  of  the  assumed  initial  conditions. 

•  The  oblique  shock  structure  predicted  within  the  passage  was  qualitatively  similar 
to  that  observed  in  the  schlieren  photographs  obtained  from  the  experiment  at  lower 
pressure  ratios. 


The  following  are  recommended  to  advance  the  experiment: 

•  Perfect  the  back  pressure  control  valve  by  sealing  and  rigidly  attaching  the  ramp 
and  replacing  or  redesigning  the  eccentric  drum. 

•  Install  pressure  taps  in  the  lower  blade  (as  designed)  and  connect  to  a  high-speed 
scanning  data  system. 

•  Manufacture  additional  windows,  or  provide  aluminum  blanks  for  tests  not  requiring 
optical  access. 

•  Provide  a  tight  shut-off  valve  in  series  with  the  present  tunnel  control  valve. 

•  Provide  a  safety  rupture  diaphragm  upstream  of  the  back  pressure  control  valve. 

The  following  are  recommended  to  advance  the  numerical  simulation: 

•  Ensure  grid  spacing  continuity,  orthogonality  along  the  boundary  of  the  grid  and 
improved  resolution  in  the  boundary  layer  region  for  use  with  the  viscous  solver. 

•  Investigate  the  "starting"  process. 

•  Run  the  Euler  solver  with  the  flow  incidence  at  1.15°  to  the  suction  surface. 
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TABLE   IA.      Run    001    Pressure   Data,    Pj/pj    =   0.765 


Port  Pressure  (psia) 

1  14.8 

2  50.0 
4  40.3 
6  17.9 
8  13.7 


TABLE  IB.   Run  001  Pressure  Data,  p2/Pi  =1.39 


Port  Pressure  (psia) 

1  14.8 

2  50.2 
4  22.4 
6  14.3 
8  19.8 


TABLE  IC.   Run  006  Pressure  Data,  p2/Pi  =1.97 


Port  Pressure 

1 

2 

4  55.8 

6  18.1 

8  35.7 
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Figure  1.   Shock-Induced  Separation  on  a  Fan  Blade  [Ref.  1] 
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Figure  2.   Vortex  Generator  Configurations  [Ref.  2] 


25 


CD 

ro 

lO 

CD 

> 

O 

C-» 

CD 

cu 

CN| 

en 

CU 

i_ 

Z3 

-  ^2 

CO 

CO 

\          a> 

CU 

\       •  — 

\      ^3 
\      f 
\     d 

_  CD 

CO 

cz 

—3 

\  DD 

X 

o 

o 

CJ3 

> 

CU 

CD 

\ 

> 

\               \ 

-  lO 

0£ 

\               \ 

c: 

\                 \ 

o 

\                  \ 

o 

X                \ 

-  CD 

o 

I— 

CU 

-«— ' 

c 

) 

_  UD 
i 

CD 

a 

1 
3          f*. 

1          1 

1              1              1 
to       ^-       r*- 

1         1 

r 

~         C 

i 

D 

c 

6        CD 

CD         CD 

CD        CD         C 

>        CD         CD         CD         CD 

Hd/d 

Figure  3. 


Typical  Surface  Pressure  Distribution 
in  the  Interaction  Region 


26 


a  =15?45° 

P  =0? 30 "BO1;  90' 

Do        (toP view)  nn       1  " 

Jet  orifices    -vl/  32 

V/p 

Flow*  -r-o-1-11 

X       /  Flow 

^VGJ  + 

llJ^jiflUIIUltlilllllllll 

(side  view) 


Figure  4.   Vortex  Generator  Jets  Configuration  [Ref.  2] 
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Figure  5.   Cascade  Geometry 
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Figure  6.   Schematic  of  Transonic  Cascade  Wind  Tunnel 
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Figure  7 .   RHS  of  Transonic  Cascade  Wind  Tunnel 
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Figure  8.   LHS  of  Transonic  Cascade  Wind  Tunnel 
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Figure  9.   Test  Section  Interior  (w/o  blades) 


Figure  10.   Schematic  of  Back  Pressure  Valve 
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Figure  11.   Back  Pressure  Valve 
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Figure  12.   Control  Valve 


Figure  13.   Control  Panel 
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Figure  14.   Schematic  of  the  Data  Acquisition  System 
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Figure  15.   Data  Acquisition  System 
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Figure  16.   Schematic  of  Optical  System 
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Figure  17.   LHS  of  Optical  System 
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Figure  18.   Schematic  of  Test  Section 
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Figure  19.   Test  Section  Inner  Wall  and  Center  Blade 
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Figure  20.   Test  Section  Upper  Blade 
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Figure  21.   Test  Section  Lower  Blade 
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Figure  22.   Test  Section  Inner  Plate  and  Model  Blade 
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Figure  23.   Test  Section  Side  Piece 
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Figure    24.       Run    001    Schlieren   -   p2/p!    =    0.765 


Figure   25.      Run   001    Schlieren   -  p2/p!   =1.39 
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Figure  26.   Channel  C-grid 
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Figure  27.   Cascade  C-grid 
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Figure  28.   Detailed  View  of  Cascade  C-grid 
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Figure  29.   Detail  of  Leading  Edge  of  Cascade  C-grid 
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Figure  30.   Detail  of  Trailing  Edge  of  Cascade  C-grid 
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Figure  31 .   Multiple  Plot  Cascade  C-grid 
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Figure  32.   Detailed  View  of  Boundary  Area 
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Figure  33.   Channel  Flow  Mach  Number  Contours 
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Figure  34.   Cascade  Flow  Mach  Number  Contours 
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Figure  35.   Cascade  Flow  Pressure  Contours 
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Cascade  Blade  Surface  Cp  Distribution 
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APPENDIX  A 
MACHINE  DRAWINGS  OF  TEST  SECTION  COMPONENTS 

Drawing  No.  Title 

T001  Test  Section  Side  Piece 

T002  Test  Section  Inner  Frame 

T003  Test  Section  Inner  Frame  and  Blades 

T004  Detailed  View  of  Side  Piece  Slot 

T005  Side  Piece  Window  Retainer  Frame 

T006  Side  Piece  Window 

T007  Sectioned  View  of  Window  Installation 

T008  Detailed  View  of  Window  Slot 

T009  Test  Section  Top  Blade 

TO  10  Test  Section  Middle  Blade 

TO  11  Test  Section  Lower  Blade 

TO  12  Pressure  Tap  Locations  on  Lower  Blade 
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Figure  Al .   Test  Section  Side  Piece 
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Figure  A2 .   Test  Section  Inner  Frame 
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Figure  A3.   Test  Section  Inner  Frame  and  Blades 
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Figure  A4 .   Detailed  View  of  Side  Piece  Slot 
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Figure  A5 .   Side  Piece  Window  Retainer  Frame 


62 


13] 


J,_ 


3 

o 

7 

i 


-3 


fN- 


s 


-fc»i 


/ 
// 


/ 


y_  >L 


r 


( + 


+* 


*  >: 


■l 

■i 

I 

,1 

i 

,_ 

2 

1,1 

e_ 

± 

3 

U) 

a 

>b 

D 

1 

.< 

in 

~~ 

I 

i> 

0 

J 

<i 

l/t 

V 

'I 

o 

u\ 

1 

z 

M 

i 

_1 

il 

<f 

0 

<i- 

■3 

'J 

7 

'J 

7 

"2 

-1 

4 

-1 

t 

2 

3 

3 

ci 

OJ 

>ll 

1 

t 

0 

£ 

-2 

V 


Figure  A6.   Side  Piece  Window 
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Figure  A8 .   Detailed  View  of  Window  Slot 
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Figure  A9 .   Test  Section  Top  Blade 
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Figure  A10.   Test  Section  Middle  Blade 
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Figure  All.   Test  Section  Lower  Blade 
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Figure  A12.   Pressure  Tap  Locations  on  Lower  Blade 
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APPENDIX  B 
DATA  ACQUISITION  PROGRAM 


10 

REM *****TXC001***** R.P.SHREEVE 3  14  *1 

20 

REM 

30 

REM   THIS  PROGRAM  I?  FOP  DATA  ACQUISITION 

40 

REM   FROM  R  SCAHI VALVE  ON  THE  TRANSONIC  CASHflDE 

50 

REM 

60 

REM   VARIABLES: 

70 

REM 

80 

REM 

90 

REM 

100 

REM 

110 

REM 

120 

REM 

130 

REM 

140 

REM 

150 

REM 

1  KM 

REM 

170 

REM   DOTE:   MRR  14-  1990 

lRfl 

REM 

190 

H I M  R[  5-10].  Ii[  5  -  1 0  ]  -  F'[  1 0  ]  •  C*[  26  3 >  R*[  f.  J 

200 

MRT  A  =  ZEP 

210 

MRT  D=ZER 

220 

MRT  P=ZER 

230 

PI  1  ]=1 

231 

Ft  2  3  =  2 

232 

P[?]=3 

234 

P[4]=4 

235 

PC5]=5 

237 

Pt 6 3=6 

233 

R[  7  ]=7 

23Q 

PC  8 3=3 

i   240 

DISP  "ENTER  ATMOSPHERIC  PRESS  <IH  HGV'5 

i   250 

INPUT  P8 

j   26" 

DISP  "ENTER  NO.  POINTS  TO  TAKE  <NAX=5>"; 

■       '270 

INPUT  SI 

£.1  -J 

DISP  "PRESS  <COHT>  WHEN  READY"! 

]   276 

STOP 

I   230 

FORMAT  B 

2?0 

FORMAT  2B 

!   300 

FORMAT  4B 

!   310 

FORMAT  F3.0 

320 

REM 

330 

REM   INITIALIZE  BUS  AND  ASSIGN  DATA  ELEMENTS 

340 

REM 

350 

IIP  I  TE  (13-  30O  >  256  -  20 .  768 ,  5 1 2  5 

360 

CMD  "?D#"» "F1R7M3A0H0T3" 

370 

REM 

330 

REM 

330 

REM 

400 

V=l 

410 

C1=0 

420 

S=l 

430 

REM 

440 

REM  BEGIN  DATA  SAMPLING 

450 

REM 

460 

FOR  K=l  TO  SI 

470 

FOR  J=l  TO  8 

Figure  Bl.   Data  Acquisition  Program  TXC001 

70 


47b  LIU1  'VI'' 

4  76  WRITE  (13- 280) V? 
480  n=Pt  J] 

490  GO SUB  1450 

Son  T 1=0 

510  REM 

520  REM  PRESSURE  DATA 

5?0  REN 

54  0  i"  Ml"  "  TD'  " 

550  MR1  IE  •:  1 3 1- 280) V 5 

568  FOR  L=l  TO  5 

5  70  i  Mli  "?D!  " 

500  HP  I IE  (13- 310)V+9 
500  HMD  "?D#"» "T3" 
€00  CMD  "•"Ct" 
610  EM1ER  '■  13  j*  '  I •  1 
620  T1=T1+DJ 

6?o  next  l 

640  Al=Tl-'5 

650  l"l  I  <  J  3=01 

651  HI  I  <  10  ]=P8 

660  me::t  j 

k70  GOTO  770 

680  REM RERU  RTMOS  TEMPERATURE  OH  SCANNER  CHANNEL  0 

600  REM 

700  CMD  "?H("  -  "C" 

710  OUTPUT  ( 13- 310  'CI 

720  CMD  "■"Dtt"-  "TO" 

730  CUD  "•"'CI" 

740  ENTER  <  13j*)DH:  -  7] 

750  CMD  "?D("> "C" 

760  REM 

770  REM 

780  PRINT 

700  PRINT  "DfiTR  SRMPLES  FOR  POINT"! 

800  PRINT 

810  FOR  1=1  TO  10 

820  WRITE  (15-830)1  -  Dll  -  I  3 

8  JO  FORMAT  F3. 0 1 SXj F 13. 6 

840  ne::t  i 

850  disp  "is  data  valid?  1=yes  0=retrke  "5 

860  INPUT  R2 

870  IF  R2=l  THEN  890 

830  GOTO  470 

890  ne::t  k 

900  REM 

910  REM  DEFINITION  OF  DATA  ARRAY 

920  REM 

940  REM  D(*)l)  =  OPEN  (RTMOS  PRESS.) 

950  REM  D(*>2)  =  CALIBRATION  PRESS. (SNITCHED) 

960  REM  IK*- 3)  =  OPEN 

970  REM  D(*>4)  =  PLENUM  PRESSURE 

930  REM  IK*- 5)  =  OPEN 

990  REM  D(*»6)  =  UPSTREAM  STATIC  PRESSURE 

lOOO  REM  IK*.  7)  =  OPEN 

1010  REM  D(*i8>  =  DOWNSTREAM  STATIC  PRESSURE 

1020  REM  IK*.  9)  =  OPEN 

1O30  REM  IK** 10)  =  ATMOSPHERIC  PRESSURE 

\  n  -i  n   p  EM**  $■■$■**  *  ■$■  %■$■$  $  $■  #  ■&#  ■&■$  ■$■£■  ■&  -s  ■$  ■$  s-  ■&  •$■  ■$  ■$■  $■  #  *  *  *  * 

1050    REM 

1O60  REM  DATA  CONVERSION  ?•  STORAGE 

107©  REM 

1O30  FOR  (=1  TO  SI 

1 090  At  K  -  1  3=Dt  K  -13*1  OOOOO 

1100    AIK.2  3  =  DIK.2  3*1OOOO0 

1110    RE  K > 3 3=Dt  K , 3 3*1 OOOOO 

1 1 20    Al  K -  4  3  =  Dt  K  -43*1  OOOOO 


Figure   Bl     (cont) .      Data   Acquisition   Program  TXC001 
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1  1  39 

Ml  \   •  5  ]=D[  I.51-  10RO0O 

]  140 

ML  1  . 6  J  —  1 1 r 1  • 6  ]* 1 00000 

1  150 

ML  1  .«7J-Li[t  t  7  ]*1 000*34.  ?2??+?.2    M4^ 

1  160 

fit  K-8]=D[Kj8] 

1 1  70 

MCI  •  P  D  ^  Lt  C  J  >9  ] 

1  180 

ML  1  -  10]  =  D[  |  .  10] 

1  190 

mf-:::  l  i 

1  1  95 

disp  "open  data  file  nun  com.  1200"; 

1  196 

SI  OP 

1  2  0  0 

DISr  "ENTER  DATA  FILE  NAME"! 

1210 

input  Mr 

1220 

FILES  * 

1  2  30 

ASSIGN  Mr.  1.1  9 

1240 

MAT  PRINT  «  i:n 

1250 

PPIMI  "DATA  stored  IN  FII  E  "Mr 

1260 

PRINI 

1270 

PEN 

14  35 

PRINT 

14  36 

PRINI  "ATMOSPHERIC  PRESSURE  =  "PR 

14  40 

END 

1450 

REN***************************************** 

1460 

REM   SUBROUTINE  "POSIT" 

14  70 

GO SUE  1650 

1480 

D=fl-P 

1 490 

HMD  "?D!  ■• 

150O 

IF  D  0  I  MEN  1580 

1510 

if  n:  n  then  1580 

1520 

RETURN 

15  30 

REN     HUME  S-V 

1540 

WRITE  ',  13.310  >V+4 

1550 

WRITE  <13j#)"C" 

1  560 

WAIT  4  00 

1570 

GOTO  14  70 

1530 

REM      ADVANCE  -  V 

1  6  0  0 

WRITE  U  3.  310  •'•/- 1 

1  6  i  0 

WRI IE  < 13j*)"C" 

1620 

WRIT  Zi) 

1640 

GOTO  1470 

1650 

PEN     READ  8  V  ADDRESS 

■   1660 

CNH  "?G*" 

1670 

P0=RBYTE13 

1680 

L=BIflHD<P0> 15) 

1690 

T=ROT' P0.4  > 

1700 

M=BIRHD<T>7> 

1710 

P=10*N+L 

1720 

WRITE  <13j2S0>256>955 

1730 

RETURN 

Figure  Bl  (cont) .   Data  Acquisition  Program  TXC001 
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APPENDIX  C 
PROGRAM  TO  DRAW  A  WEDGE-ARC  CASCADE 


10  REM *****SSCB03****s R. P. SNREEVE 

10  19  in 

20  REN PROGRflM  TO  CRLC.  ?■  I'RFil J  (h  UEDGE-RRC  CRS 

;nDE  BLFIDE 

30  REM INTO  n  CRSCRBE 

35  H EG 

40  URTR  6. 0. 015. 3. 5.2. 85.-51. 34. 2- 30 

50  reru  CO •  RO «  ft  1 .  11 1  >  GO  >  SO .  G9 

51  R0=-G8 

55  I! 1=R0*< 1 -C OSRO-S I HhO > 

56  D2=R0*<  1  -COSRO+S I NR0 ) 

57  S1=C0*S1NR0  SO 

53  S2  =  CP  SO 

60  REIT CflLC.  BLHUE  GEOM. 

ro  :'o=po 

30  V0=R0 

i6o  xi=x9-R0*sinni 

1 10  Y1=R0+R0*COSR1 

120  ;:2=ro 

130  V2=0 

140  E0=<R0  inn-  Ml  2)}-Y.2 

150  R8=C0-R0-U1 

160  R3=R8*R8 

170  B8=R0-<Ul+E0>*TRHfil 

130  B9=B8*B8 

130  Rl  =  <R9+B9-R0*R0)/'<2*<:R8*SIHRl-B8*COSni-R0>> 

200  X3=W1 

210  Y3=<Ul+E0>*TRHfll 

220  K4=X3+R1*SINR1 

230  V4=Y3-R1*C0SR1 

240  ''.7  =  C0-R0 

250  Y7=R0 

260  P1=RTM',  <M7-;M>,'<Y7-Y4  >> 

270  X5=K7+R0*SINP1 

230  Y5=Y7+R0*COSP1 

230  ';'.h  =  VJ 

300  Y6=0 

310  REM FRUIT  RESULTS 

320  PR  HIT 

330  PRINT 

340  PRINT  "WEDGE-ARC  CRCRDE  BLRDE  " 

350  PRINT 

360  PRINT  "CHORD="C0. "WEDGE  LENG1M="W1 

370  PRINT  "WEDGE  RNGLE="R1  "DEG"  »  "L.  E.  RRDIUS="R0 

330  PRINT  "SURFRCE  RRDIUS="R1 i "PHI="P1"DEG. " 

330  PRINT  "X0="}{0 

400  PRINT  "Xl  =  ">Ui "Y1="Y1 

410  PRINT  "!'.2="X.2."Y2="Y2 

Figure  Cl.   Cascade  Geometry  Program  SSCB03 

73 


420  r f  i m i  ";•::-"::?■ ' 

Y3= " • V 

~; 

4?0  PRIHI  "}{4  =  "!Mj  ' 

Y4="Y4 

440  PRINT  "}<5="K5.  ' 

Y5  = "  Y5 

47.0  PRINT  "}{6="M6»  ' 

Y6-"Y6 

460  print  "X7="H7»  • 

Y7= " Y 7 

4  70  PRIUT 

471  PPIHT  "cnscnuE 

GEOMETR' 

472  PPIHT  "SETTING 

FlNGLE= 

"fl0" 

1' EG.  " 

.  "BLHl't 

475  PRINT  TOP  SOLID ITY=" 

so. " 

BLADE 

spncn 

400  PRINT 

4  90  STOP 

1000  REM LiRlTUBLnUE 

1010  SCflLE  -It  14 » -J 

•*  *  x 

1020  LtlBEL  <*>  l.Sti 

1  •  0 . 10 

15) 

1030  TOP  T=-?0  TO  ' 

?0-fll   S 

TEP 

(180- 

A  l  ) 

30 

1040  J!9=R0*(  1  -COS  1 

> 

1041  Y?=R0*<l+SIin 

) 

1042  G 0 S U B  5 0 0 0 

1050  NEXT  T 

1060  ;:?=:;?• 

1061  Y9=Y3 

1062  GO SUB  5000 

1070  TOP.  T=-ni  TO 

PI  STEF 

<  n  i 

f  p  n , 

250 

ioso  ::?=:'4+ri*simt 

1001  Y9aY4+Rl*C0ST 

1002  GO SUB  50OO 

1090  NEXT  T 

1 lOO  FOR  T=-90+Pl 

1  0  90  ': 

•  TEP 

'180- 

-PI) 

'30 

1110  :V?=O;7  +  R0*COST 

1111  Y9=Y7-R0*SINT 

1112  GO SUB  5OO0 

1120  NEXT  T 

1 130  ;:?=po 

1131  Y9=0 

1 132  GOSUB  50O0 

1  140  PEN 

1150  I" ISP  "  INPUT  2 

POP  BLf) 

2.  3 

FOP 

■~.\ 

1151  IMPUT  Z0 

1 153  O1=(Z0-1 )*S2* 

C0SG9 

1154  02**<Z0-1>*S2* 

SING9 

1155  OFFSET  01. 02 

1156  GOTO  103O 

1160  END 

5000  REM SUBPL 

UTINE 

TO  TRANSLATE 

FlUD 

5010  X=<X9-D1)*C0S 

G0+<.  Y9 

-D2> 

*SING 

n 

5O20  Y=<Y?-D2)*C0S 

G0-CX9 

-Dl> 

*SING 

O 

5025  PLOT  X,Y 

5O30  RETURN 

SPACEOi; 


rotate  axes 


Figure  CI  (cont) .   Cascade  Geometry  Program  SSCB03 
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APPENDIX  D 

DESCRIPTION  OF  THE  NUMERICAL  METHOD 
A.      REYNOLDS-AVERAGED  NAVIER-STOKES  EQUATIONS 

The  Reynold 's-averaged  Navier-Stokes  equations  are  presented  first.  The  unsteady, 
compressible  Navier-Stokes  equations  are  developed  in  Ref.  9.  A  complete  derivation  of 
the  Reynold's  averaging  process  can  be  found  in  Ref.  10,  which  is  the  source  of  the 
following  description. 

The  time-dependent,  compressible  Navier-Stokes  equations  in  tensor  notation  are 
as  follows: 

Pt  +    (P"j>j  =  0 

(pui)t  +    (pUjUj   +  Oyjj  =  0  i   =  1,2,3  (l) 

(p£)t  +    {pEUj  +  uiaij  +  qj)  j  =  0 

where  subscript  notation  has  been  used  to  denote  partial  derivatives,  eg.  pt  =  ijP   and  the 

summation  convention  has  been  used  for  repeated  indices.     The  stress  tensor  for  a 
Newtonian  fluid  is  given  as 

Here  the  Stokes  hypothesis  3X  +  2  =0  has  been  used.  The  heat  flux  vector  is  obtained 
from  the  Fourier  law  of  heat  conduction  as 
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Qj   =    -KTj    =    -JLhj  (3) 


To  derive  the  Reynold's- averaged  Navier-Stokes  equations,  the  averaging  operation 
is  defined  as  follows: 


pixjt  t)    =   -±-fp(xjrs)ds 


t+T 

(4) 


2T 

C-T 


where  2T  is  the  time  averaging  interval,  which  is  assumed  to  be  large  enough  compared 
to  the  energy  containing  turbulent  time  scales  but  small  when  compared  to  the  time  scale 
of  the  average  motion.  It  is  assumed  that  there  is  a  range  of  values  of  T  where  the  mean 
density  is  independent  of  T. 

The  fluctuating  density  is  defined  as  the  difference  between  the  density  and  its  mean 
value. 

pUp.p  (5) 

The  averages  and  fluctuating  values  for  the  other  variables  are  defined  similarly.  It  is 
assumed  that  the  average  of  the  fluctuations  equals  zero  and  the  average  of  the  average 
still  is  the  average. 

The  time-averaged  Navier-Stokes  equations  are  obtained  by  averaging  Eq.  (1).  The 
resulting  equations  involve  terms  similar  to  the  following: 
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pUtUj   =    p    Ut    Uj    +    p    U'  iU'  j    +    U±    P'U'j   +    Uj    p'u'j,    +    p'u'jU'j  (6) 

For  incompressible  flows,  where  the  average  of  the  density  fluctuation  equals  zero, 
the  last  three  terms  of  Eq.  (5)  are  dropped.  The  compressible  equations  contain  additional 
terms,  such  as  uj1"  p'  u/  ,  therefore,  an  alternative  method  for  averaging  velocity  and 
energy  terms  is  followed.  This  method  is  called  mass-averaging  (Favre-averaging)  with 
the  following  definitions: 

ui  =  pu~/p  K  =  p72/p  <7> 

u' d  =  ud  -  ui  h%  =  h  -  K  (8) 

It  is  important  to  note  that  the  averages  of  the  fluctuating  quantities  for  Favre-averaging 
is  not  zero. 

uTT~1  =  -p'u'j/p  ~ErT  =  -p'A'/p  <9) 

However,  the  mass-weighted  variables  are  zero. 


pu«  'j  =  0  ph<  •    =  0  dO) 

By  introducing  the  mean  and  fluctuating  quantities  into  Eq.  (1)  and  averaging,  the 
Reynold's-averaged  Navier-Stokes  equations  are  obtained. 

p^  +    (pUj)^  =  0 


(P"i)  t  +    (P  uiuj  +  °ij)j  =  ° 


(ID 
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ipE)  t  +    {pEQj  +   GjOij  +  Qj)j  =   0  i  =  l,2,3 

0iJ.  =  0,/  +  o,/,  gj  =  g/  +  g/ 

The  following  are  identified  by  name: 

Reynold's  Stress  Tensor:       a1:fT  =  pu±'  •  u^ •  ' 


o 

Viscous  Stress  Tensor:  a±jM  -  -pfijj  -  M"^  +  Qj.i  -  -^-&ijc*./c) 


(12) 


Reynold's  Heat  Flux  Vector:  g^  -  ph%  Hip  ' 

Molecular  Heat  Flux  Vector:  g^M  =  --^-R^ 

Total  Energy:  E  =  e  +  k  +  ukuk/2 

Turbulent  KE:  k  =  pu^'  ■  u*'  '/2p 

It  is  assumed  that  (1  is  independent  of  time,  but  depends  on  temperature  variations. 

The  Reynold's  stresses  are  related  to  the  mean  flow  quantities  through  turbulence 
modeling.  The  turbulence  models  used  in  the  present  work  are  discussed  in  the  next 
section.  When  the  Reynold's  stresses  and  heat  fluxes  arc  related  algebraically  to  the 
mean  flow  quantities,  the  resulting  models  are  called  algebraic  turbulence  models.  The 
algebraic  turbulence  models  are  either  mixing  length-type  models  or  eddy  viscosity-type 
models.  In  the  present  study  the  widely  used  Baldwin-Lomax  algebraic  eddy-viscosity 
model  was  used. 
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The  eddy  viscosity  models  are  zero-equation  models,  and  the  eddy  viscosity  is 
modeled  algebraically.  More  complex  models  exist  that  take  into  account  the  transport 
of  certain  turbulence  quantities  such  as  kinetic  energy  or  shear  stress.  These  models  are 
derived  from  simplified  forms  of  the  Navier-Stokes  equations  for  mean  and  fluctuating 
quantities  [Ref.  11]. 

A  simple  one-equation  model  that  is  specially  designed  to  predict  shock-boundary 
layer  separated  flows,  the  Johnson-King  model  [Ref.  12]  was  used  in  the  present 
investigation.  This  model  was  considered  to  be  most  suitable  for  the  computation  of  the 
nonequilibrium  flow  in  the  recovery  region  after  the  shock-induced  separation,  based  on 
the  experimental  data  given  in  Ref.  12.  Figure  Dl  shows  a  comparison  of  predictions  for 
the  separated  flow  region  after  a  transonic  shock  using  four  turbulence  models.  The 
Johnson-King  and  RNG  models  exhibit  the  best  results.  The  Baldwin-Lomax  model  and 
the  Johnson-King  model  are  briefly  described  in  the  next  section. 

B.      TURBULENCE  MODELING 

1.       Baldwin-Lomax  Model 

The  Baldwin-Lomax  model  is  a  two  layer  algebraic  model  which  does  not  require 
the  finding  of  the  boundary  layer  quantities,  as  does  the  Cebeci-Smith  model.  Instead, 
the  outer  and  inner  eddy  viscosities  are  scaled  by  the  vorticity.  A  complete  formulation 
of  the  model  is  given  in  Ref.  13. 


79 


2.       Johnson-King  Model 

The  shock- separated  region  is  computed  with  the  Johnson-  King  turbulence  model. 
A  detailed  description  of  the  development  of  the  model  is  given  in  Ref  12.  The  basic 
formulation  of  the  model  necessary  for  the  numerical  implementation  is  presented  briefly. 

The  Johnson-King  turbulence  model  was  developed  for  the  purpose  of  improving 
the  predictions  of  pressure  driven  separated  flows  and  transonic  shock-induced  separated 
flows.  The  Cebeci-Smith  model  [Ref.  10]  was  proven  fairly  accurate,  even  though  it  was 
based  on  the  invalid  assumption  that  the  turbulent  shear  stress  depends  only  on  the  local 
properties  of  the  mean  flow.  The  Johnson-King  model  develops  a  closure  for  boundary- 
layer  flows  without  attempting  to  predict  the  production,  dissipation  and  diffusion  of  the 
turbulent  kinetic  energy  and  the  dissipation  rate,  throughout  the  entire  viscous  layer,  as 
in  two-equation  models. 

To  account  for  convection  and  diffusion  effects  on  the  Reynolds  shear  stress 
development,  an  ODE  is  developed  to  describe  the  maximum  Reynolds  shear  stress 
development  in  the  streamwise  direction.     The  ODE  describes  the  development  of 


u'  v'  =  y^xt/p    along  the  path  of  maximum  shear  stress. 
The  eddy  viscosity  is  given  by 

vt 
vt  =  vt  [1  -  exp( '-)]  (13) 

where  vu  and  vto  describe  the  eddy  viscosity  variation  in  the  inner  and  outer  parts  of  the 
boundary  layer.   The  inner  eddy  viscosity  is  computed  as 
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D2  =  1  -  e{~z/An 
where  the  constant  A+  =  15.   The  outer  eddy  viscosity  is  given  by 

vt    =  o(x)  •[0.0168C7e8*Y]  (15) 

where  7  is  Klebanoff  s  intermittency  function 

v  =    [1   +  S.S(i)6]-1  (16) 

o 

and  o(x)  is  the  solution  of  the  ODE  mentioned  earlier. 

To  complete  the  formulation,  an  equation  is  needed  for  w  v  a.  This  equation  is 
developed  from  the  turbulence  kinetic  energy  equation.  The  result  is  the  following  ODE 
after  a  change  of  variables: 

dx       2^lU       gj        a.atO.7  -  (y/8)J  ]1       X  J    |J   ' 

where  C^f  and  a,  are  modeling  constants,  ua  is  the  maximum  average  mean  velocity,  g 
and  g^  are  given  as 


g  =  l>y(-u'v')„     and     g^  =  l/^-u-  V  )  m>eQ  U8> 

and  L,,,  is  the  dissipation  length  evaluated  as 

Lm  =  0.40z         zjb<  0.225  (19) 

L,,,  =  0.095        zJS  >  0.225 
The  equilibrium  shear  stress,  g^,  is  determined  from  the  following  equilibrium  eddy 
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viscosity  distribution 


Vtj.eg  =   iPxZyH-lHv*)^  <2°> 

vto,e<7=  0.0168ue6*Y 
An  implicit  Euler  method  is  used  for  the  numerical  solution  of  Eq.(17),  and  the 
maximum  shear  stress  at  each  iteration  level  is  updated  as  follows: 

vf1 
a(x)n+1  =  a(x)n — -  (2D 

Solutions  with  the  Johnson-King  turbulence  model  were  obtained  as  follows.  First  a 
convergent  solution  using  the  Baldwin-Lomax  turbulence  model  for  the  entire  flow  field 
was  obtained.  Then  the  Johnson-King  model  was  applied  to  the  upper  part  of  the  airfoil. 
To  initiate  the  solution  a(x)  in  Eq.(15)  is  set  equal  to  one,  and  it  is  allowed  to  change 
according  to  Eq.(21)  until  the  final  solution  is  obtained.  Note  that  the  Johnson-King 
model  reduces  to  the  Cebeci-Smith  model  when  o(x)  is  identically  equal  to  one. 

C.      COMPUTATIONAL  METHOD 

The  following  section  was  developed  using  notes  from  J.  Ekaterinaris. 

1.      Discontinuous  Solutions  and  Entropy  Method 

The  numerical  solution  of  the  initial  value  problem  for  hyperbolic  systems  of 
conservation  laws  consists  of  the  solution  of 
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Ut    +    fx(u)     =    0  (22) 

u{x,  0)    =   u0(x) 

where  u(x,t)  is  the  vector  of  m  unknowns  and  f(u)  is  the  flux,  a  vector  valued  function 
of  m  components.   The  matrix  form  of  Eq.(22)  is 

"t  +  Aux  «  ° 
Because  the  system  is  assumed  hyperbolic,  all  the  eigenvalues  of  the  Jacobian  matrix  are 
real  and  positive. 

To  allow  for  discontinuous  solutions  we  admit  weak  solutions  which  satisfy  Eq.(22) 
in  the  sense  of  distribution  theory  as  shown  below 

f    f  {<f>tu  +  $xf [u))dxdt  +    [$(x,0)u0(x)dx  =  0         (23) 

z=0;t=-«»  -" 

where  0  are  C~  test  functions  that  vanish  for  |x|  +  t  large.  The  hyperbolic  system  is 
satisfied  in  a  pointwise  sense  by  a  piecewise  smooth  solution  in  each  smooth  region, 
while  across  each  discontinuity  the  Rankine-Hugoniot  relation 

f(uR)    -  f(uL)    =  C(uR  -  uL)  (24) 

holds  where  C  is  the  speed  of  propagation  of  the  discontinuity,  and  the  left  and  right 
states  are  designated  uL  and  uR,  respectively. 

The  conservation  law  Eq.(22)  is  considered  to  possess  an  entropy  function  E(u) 
defined  as 
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•  E  is  a  convex  function  of  u,  Euu  >  0  convexity 

•  E  satisfies  Eufu  =  Fu  compatability 

where  F  is  the  entropy  flux  function. 

The  convexity  condition  forces  irreversible  processes  to  run  in  the  correct  direction 
and  produce  entropy.  A  consequence  of  the  compatability  condition  is  that  the  reversible 

processes  do  not  produce  entropy.    Note  that  the  mathematical  entropy  E,  and  physical 

i 
entropy  S  =  logCp/p1'),  statements  are  compatable,  and  if  one  is  true  then  the  other  is  true. 

Every  smooth  solution  of  Eq.(22)  then  satisfies 

Eu(t)    +  Fx(u)    =  0  (25) 

If  u  is  piecewise  smooth  with  discontinuities,  then  Eq.(25)  holds  pointwise  in  the  smooth 
regions;  and  across  a  discontinuity 

F(UR)    -  F(uL)    <L  C[E(UR)    -  E(uL))  (26) 

Consider  for  simplicity  that  the  numerical  approximation  to  weak  solutions  of 
Eq.(22)  are  obtained  by  the  following  explicit  scheme 

uT1  =  "i5-  £(#1/2  "  #W  (27) 

where  f"tl/2  is  a  numerical  flux  defined  as 

It  is  required  that  the  numerical  flux  is  consistent  with  the  physical  flux  in  the  following 
sense 
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f(u,u)    =  f(u) 

In  addition,  the  differencing  scheme  is  consistent  with  the  entropy  condition  if  the 
inequality  of 


£j+1  *  E?  -  j*{Ff.l/2  +  F?_1/2)  (28) 


is  satisfied  where 


E?1  =  EiuT1) 
and 

Fi"l/2  =  F(u?,u?+1) 
Assume  that  the  difference  scheme  in  Eq.(27)  is  consistent  with  the  hyperbolic 
conservation  law  Eq.(22)  and  the  entropy  condition  Eq.(26).    Let  u°  be  a  solution  of 
Eq.(27)  for  initial  values  u°  =  fl(j'A) .  Suppose  that  for  lim  Ax  — »  0  (Ax  =  x/i^)  and  At 
->  0  (At  =  T/N) 

limAx-o  u"  =  u[x,  t) 

At-0 

the  limit  exists  in  the  sense  of  L,  norm  convergence.  Then  the  limit  satisfies  the  weak 
form  of  the  conservation  law  Eq.(23)  and  the  weak  form  of  the  entropy  condition  given 
by 
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-ff(4>tE  +  bxJEDdxdt  -   [$(x,0)E{u0{x)  )dx  <>  0 

o—  -- 

where  4>  is  a  nonnegative  smooth  test  function  <}>(x,t)  possessing  a  compact  support. 
Assume  that  the  initial  data  corresponds  to  some  reference  state  u*  for  large  |x| 

u0(x)    =  u*  \x\  >  M 

then 

u"  -  u*      for      Axi  >  M  +  nAx 

Without  altering  the  convexity  of  the  entropy  function  an  inhomogeneous  function  can 
be  added  to  obtain 

E(u')    =  0        EUi(u')    =  0 
By  summing  Eq.(28)  with  respect  to  i  to  obtain 

£    Ef1  s.  £    Ef 

i  i 

shows  that  entropy  is  a  decreasing  function  of  time,  and 

i  i 

indicates  that  the  scheme  is  stable. 

This  condition  is  not  strong  enough  to  prove  pointwise  boundness  of  the  solutions 
or  the  existence  of  convergent  subsequences.  Solutions  of  equations  with  compressible 
flow  must  ensure  that  the  difference  scheme  used  keeps  the  variables  within  their  physical 
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range  (density  and  pressure  are  always  positive  quantities).  The  above  development  is 
true  for  any  number  of  space  variables.  In  addition,  they  hold  for  multidimensional 
schemes  consisting  of  one-dimensional  factional  steps  provided  that  each  individual  one- 
dimensional  step  satisfies  the  entropy  inequality. 

2.      Upwinding  Schemes  for  Linear  Hyperbolic  PDE's 

A  linear  constant  coefficient  hyperbolic  PDE  has  the  form 

|H  +  a|H  =  o  (29) 

dt  ox 

The  numerical  approximation  of  the  solution  at  x  =  iAx,  t  =  nAt  is  u"  where  Ax  is  the 
spatial  mesh  size  and  At  is  the  time  step.  A  simple  first-order  accurate  upwinding  scheme 
that  takes  into  account  the  direction  of  wave  propogation,  according  to  the  characteristic 
theory  of  hyperbolic  PDE's,  which  may  be  used  for  the  solution  of  Eq.(29)  is  as  follows 

[u?  -  u^  a  >  0 

where  X  =  At/Ax.  This  numerical  scheme  can  be  written  in  the  following  single  equation 
form 

if/1  =  u?  -  |a[Uin+1  -  U;.J   +  jlaKuJU  -  2un±  +  uU)    Pia) 
or 

uT  =  uj3  -  i.[a+{u?  -  uili)   +  a'(u?+1  -  u]2)  ]  Oiw 

where  a+  =  1/2  (a  +  |a|),  a  =1/2  (a  -  |a|). 
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Higher  order  upwind  schemes  may  be  obtained  by  increasing  the  order  of  accuracy  of  the 
stencils  of  Eq.(30)  in  the  upwind  direction. 

The  preceding  ideas  may  be  generalized  for  a  system  of  hyperbolic  equations  with 
constant  coefficients  of  the  following  form 

I?  +  a|?  =  o  02, 

dt  ox 

where  U  is  a  vector  of  m  unknowns  and  A  is  a  m  x  m  constant  coefficient  matrix  having 
real  eigenvalues.  Let  X  be  the  eigenvalue  transformation  diagonalizing  A  as  X  'AX  =  A, 
where  A  is  the  diagonal  eigenvalue  matrix.  Then  Eq.(32)  can  be  rearranged  in  diagonal 
form  as 

|?  ♦  A#  =  0  (33, 

dt  ox 

where  W  =  X'U. 

Application  of  the  numerical  schemes  shown  in  Eq.(31)  for  the  vector  equation  in 
Eq.(33)  yields 

WT1  =  W?  -  \MWli  -  Wlx)   +  ||A|{»ff+1  -  2W?  -  W?.,}    (34) 
This  scheme  can  be  written  in  terms  of  the  original  variables  as 

C/T1  =  US  -  \MUlx  -  U2_x)   +  \\A\{U^  -  2Un±   -  vW\      ^ 
where  |A|  =  XAX'1.   Similar  to  the  one  equation  case  the  new  variables  A+  and  A"  may 
be  defined  as  A+  =  1/2(A  +  |A|)  and  A"  =  1/2(A  -  |A|).   Equation  (35)  now  takes  the 
form 
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U"/1  =   Un±    -  X{Fi+1/2  -  F±_l/2} 


(36; 


where 


Fi+1/2  =   \[A{Ui+1  +   17,)    -    |A|(C7i+1  -    U±)} 

Fi-i/2  =  \IMU±+   U±_x)    -    \A\{U±  -   U±_x)] 
The  term  P  is  referred  to  as  the  numerical  flux. 

3.      Conservative  Schemes  and  Shock-Capturing  Theory 

The  theory  of  modern  shock-capturing  numerical  schemes  for  non-linear  hyperbolic 
conservation  laws  such  as  the  Euler  equations  relies  on  the  basic  first-order  upwinding 
methods  discussed  in  the  previous  section.  High  order  shock-capturing  methods  suitable 
for  gas  dynamic  hyperbolic  equations  and  nonconstant  coefficient  cases  use  non-linear 
differencing  schemes.  This  significantly  affects  the  stability  of  these  methods  because  the 
stability  checking  can  only  occur  locally  and  for  the  linearized  versions  of  the  non-linear 
equations.  This  local,  linearized  stability  is  neither  necessary  nor  sufficient  for  non-linear 
problems  such  as  shock  discontinuities.  Traditional  methods  of  removing  non-linear 
instability  involve  the  introduction  of  linear  or  non-linear  numerical  dissipation  (artificial 
viscosity  or  smoothing)  into  the  difference  schemes.  This  approach  alone  does  not 
guarantee  convergence  to  a  physically  correct  solution,  and  higher  order  accurate  upwind 
schemes  have  been  developed  recently  to  overcome  the  deficiencies  of  the  artificial 
dissipation  approach. 
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The  limit  solution  of  any  finite-difference  scheme  in  a  conservation  form  consistent 
with  the  conservation  laws  satisfies  the  jump  conditions  across  a  discontinuity 
automatically  [Ref.  14].  Weak  solutions  (solutions  with  shocks  and  contact 
discontinuities)  of  hyperbolic  laws  are  not  uniquely  determined  by  their  initial  values 
except  if  the  entropy  condition,  given  in  the  previous  section,  is  satisfied  in  order  to  pick 
out  the  physically  relevant  solution.  Finite-difference  approximations,  however,  always 
converge  to  a  physically  relevant  solution  when  the  numerical  scheme  used  is  monotone. 
If  the  scheme  is  not  monotone  then  it  must  be  consistent  with  the  entropy  inequality  to 
assure  convergence.  These  monotone  schemes  possess  many  desirable  qualities  required 
for  the  calculation  of  discontinuous  solutions. 

a.      Monotone  Schemes 

Consider  the  scalar  hyperbolic  conservation  law 

Ut   +    aUx   =    0  (37) 

where  a(u)  =  J*£  is  a  characteristic  speed.  A  three-point  explicit  scheme  in  conservation 
form  for  the  numerical  integration  is 

U^   =    U?  -   AU&1/2    -   h?_1/2)  (38) 

where  X  =  Ax/At  and  the  numerical  flux  function  h",U2  ■  h(u",  u"tl)  is  consistent  with  the 
conservation  law  in  the  sense  that  Mu^Uj)  =  f(Uj). 

The  numerical  scheme  in  Eq.(38)  is  of  the  following  generalized  functional  form 
and  it  is  monotone  if  g  is  a  monotonically  increasing  function  of  each  of  its  arguments. 
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UT1    =    9{UU,U?,UU  (39) 

Monotone  schemes  have  the  attractive  property  of  producing  smooth  transitions  near 
discontinuities,  but  they  are  only  first  order  accurate.  Unfortunately,  not  all  first-order 
upwind  schemes  are  monotone.  In  addition,  monotone,  first-order  upwind  schemes  cannot 
produce  solutions  for  the  flow  field  with  the  required  accuracy  using  reasonable  grid 
spacing.  Still,  higher-order  shock-capturing  schemes  required  for  practical  computations 
require  some  linear  or  non-linear  numerical  dissipation.  Unfortunately,  these  higher  order 
schemes  suffer  from  the  following  deficiencies:  they  produce  spurious  oscillations  when 
the  solution  is  not  smooth,  they  may  develop  instabilities  at  discontinuities  and  they  may 
select  nonphysical  solutions. 

Two  classes  of  modem  shock-capturing  schemes,  the  TVD  and  the  ENO  schemes, 
are  more  appropriate  for  the  computation  of  weak  solutions,  especially  when  these 
schemes  are  consistent  with  the  entropy  inequality  and  are  second-order  or  higher  in 
smooth  regions.  The  main  distinction  between  ENO  and  TVD  schemes  is  that  the  first 
can  retain  higher  than  first-order  terms  at  extrema  while  the  TVD  schemes  reduce  to 
spatially  first-order  in  these  regions.  However,  they  are  more  efficient  in  terms  or 
operations  count.   The  next  section  lays  out  the  development  of  TVD  schemes. 
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b.      TVD  Schemes 

Consider  a  family  of  five  point  difference  schemes  of  the  form 

uT1  *  xe(hZl/2  -  hT-l/2)  =  u?  -  xa  -  e)  (hxn.1/2  -  h?.1/2)        (4o 

where 

hi+1/2  =  h(uD-lrUi,  ui+1,ui+2) 
and  the  scheme  can  be  explicit  or  implicit  depending  on  the  value  of  the  parameter  8.  If 
9  =  0,  the  scheme  is  explicit;  if  6  =  1/2,  the  scheme  represents  the  trapezoid  rule;  and  if 
0  =  1,  the  scheme  is  a  backward  Euler  method. 
Defining 

n+l 


^♦i/2  ■  (i  -  o)^i/2  +  e/2j;i/2  (4D 


Eq.(40)  becomes 


UT1    =    Uj1   -    M/2i+l/2    -   h±_U2)  (42) 


which  is  of  the  form 


L^U**1  =  L2uD 
with  Lj  and  Lj  defined  as  the  following  difference  operators 

L2u  =  Uj,  +  Ml  -  0)  (h1+1/2  -  h±_uJ 
The  total  variation  of  the  mesh  un  is  defined  as 
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TV(un)    =  S|uin+1  -  u*\  =  E|Ai+1/2un|  (43) 

The  numerical  scheme  is  TVD  if  the  following  inequality  holds 

7V(un+1)    <:  TV{un) 

There  are  four  design  principles  used  in  the  construction  of  high-resolution  TVD 
schemes: 

•  The  flux-corrected  transport  of  Boris  and  Book  [Ref.  15].  This  is  a  two  step  hybrid 
scheme  consisting  of  a  combined  first  and  second-order  scheme.  It  computes  an 
update  form  of  a  first-order  scheme  and  then  filters  the  second-order  corrections  by 
the  use  of  flux  emitters  to  prevent  occurence  of  new  extrema. 

•  A  second-order  extension  of  Gudunov's  scheme  by  van  Leer  [Ref.  16]  and  Colella 
and  Woodward  [Ref.  17].  This  scheme  is  based  on  the  observation  that  one  can 
obtain  second-order  accuracy  in  Gudunov's  scheme  by  replacing  the  piecewise- 
constant  initial  data  of  the  Reimann  problem  with  piecewise-linear  or  parabolic 
initial  data.  The  slope  of  the  initial  data  is  chosen  so  that  no  spurious  oscillations 
can  occur. 

•  The  modified-flux  TVD  scheme  of  Harden  [Ref.  18].  This  is  a  technique  to  design 
a  second-order  accurate  TVD  scheme  by  starting  with  a  first-order  TVD  scheme  and 
applying  it  to  a  modified  flux.  The  modified  flux  is  chosen  so  that  the  scheme  is 
second-order  in  smooth  regions  and  first-order  at  extrema. 

•  The  numerical  fluctuation  approach  of  Roe  [Ref.  19]  depending  on  an  average 
function.  The  average  function  is  constructed  in  such  a  way  that  spurious 
oscillations  will  not  occur  by  the  use  of  flux  limiters. 

Most  of  the  above  methods  can  be  viewed  as  three-point  central-difference  schemes 

with  an  appropriate  numerical  dissipation  term  which  automatically  controls  the  amount 

of  numerical  dissipation.   This  control  of  the  dissipation  is  automatic  unlike  the  control 

of  the  numerical  dissipation  used  by  central-differencing  techniques  used  in  linear  theory. 

The  basic  idea  of  the  above  design  principles  is  to  construct  a  higher-order  scheme  which 
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has  properties  similar  to  a  first-order  TVD  scheme  so  that  spurious  oscillations  cannot  be 
generated.  In  order  to  achieve  that  limiting  procedure,  a  flux  limiter  is  used,  which 
imposes  constraints  on  the  gradients  of  the  dependent  variables  or  the  flux  function. 

The  derivation  of  the  higher-order  schemes  with  TVD  properties  can  be  achieved 
by  either  following  the  so-called  MUSCL  (Monotonic  Upstream  Scheme  for  Conservation 
Laws)  approach  of  van  Leer  [Ref.  16],  Harden  [Ref.  18],  Roe  [Ref.  19]  and  Osher  and 
Chakravarthy  [Ref.  20];  and  the  non-MUSCL  approaches  of  Roe  [Ref.  21]  and  Yee 
[Ref.22]. 

4.       Development  of  Flow  Simulation  Code 

The  numerical  scheme  used  for  the  simulation  was  developed  in  Ref.  23.  A 
summary  of  the  numerical  scheme  is  given  next. 

The  general  form  of  the  conservation  law  form  for  an  inertial  reference  frame  is 

-4-\QdV  +  in-FdS  =   (v-SV  (44) 

V  S  V 

The  main  motivation  of  the  conservation  law  form  is  to  capture  flow  discontinuities. 
Finite  volume  discetization  is  obtained  by  replacing  the  surface  integrals  by  sums  over 
the  faces  of  a  cell.  The  differential  form  is  obtained  by  applying  Eq.(44)  to  a  cell  in 
physical  space.  In  order  to  facilitate  the  numerical  computation  for  complex  geometries, 
the  physical  space  cell  (x,z)  is  the  image  of  the  computational  domain  cell  (£,£)  resulting 
from  the  generalized  coordinate  transformation 

x  =  x(|£,t)   ,   z-zffi&t)  (45) 
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The  transformed  coordinate  system  (^,Q  is  an  evenly  spaced  rectangular  domain 
(orthogonal)  having  area  of  unity.  This  allows  all  numerical  differentiations  to  be 
performed  with  unweighted  differencing  formulas.  In  addition,  the  conservation  law  form 
represents  the  thin-layer  approximation.  Terms  of  order  l/ReI/2  and  smaller  are  neglected, 
and  the  viscous  terms  containing  derivatives  in  the  strearnwise  direction  are  dropped. 
These  viscous  terms  are  smaller  in  magnitude  than  the  viscous  terms  with  derivatives 
normal  to  the  body  surface.  The  thin-layer  approximation  allows  a  reduction  in 
computational  time,  since  the  strearnwise  viscous  gradients  are  not  usually  resolved 
adequately  anyway.  Also,  the  thin-layer  approximation  allows  computation  of  separated 
flows  and  flows  with  normal  pressure  gradients.  The  conventional  boundary  layer 
approximation  does  not  allow  these  computations. 

The  differential  form  of  the  two-dimensional  conservation  law  in  generalized 
coordinates  is  as  follows: 


dtQ  +  d^E  +  dnF  =  Re-^S 


(46) 


where  the  transformed  vectors  are 


0=  4 


pu 
pw 

{  e ) 


E  = 


pU 
puU  +  Zj) 

PVU  +    IjD 

\{e  +  p)u  -  Ztp) 


.  f=  4 


pw 

puW  +   CxP 

pww  +  CZP 
(e  +  p)   -  CtP; 


(47; 
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1 lea' 


3*  1 


0 


V-m^u^  +    (n/3)/n2Cx 

\nn1vQ  +    {\i/3)m2Cy 

^/^m,  +    (|i/3)/n2  +    (CAu  +  CyV  +  Cz^)y 


™1     =     Cx    +Cz 


m2    =    CXUC    +    C2"C 


(48) 


and  1 1  and  W  arc  the  contravariant  velocity  components  given  by 


v  ■  "£*  +  ^z  ♦  e 


(49! 


W  ■  uCx  +  "C*  +  Ct 

In  the  above  equations,  all  geometric  dimensions  are  normalized  with  the  chord  length; 
p  is  the  density  normalized  with  the  freestream  density  p„;  u  and  w  are  the  Cartesian 
velocity  components  normalized  with  the  speed  of  sound;  e  is  the  total  energy  per  unit 
volume  normalized  with  p_a2..;  Pr  is  the  Prandtl  number;  and  k  is  the  thermal 
conductivity.  The  pressure  is  related  to  the  density  and  total  energy  through  the  equation 
Of  state  for  an  ideal  gas,  p  =  (y  -  1  )[e  -  p(u2  +  w2)/2]. 

Spatial  schemes  for  the  numerical  solution  of  the  inviscid  part  of  Eq.  (46)  are  of  the 

following  semidiscrete  form 
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dQ    _    _  Hi,k+l/2    "    H,i.k-l/2    _     Gi.k*l/2    "    Gi.k-l/2  {50) 

dt  A£  AC 


where  H  and  G  are  the  fluxes  defined  as 


1 

Hi  +  \/2.k   ~    ~5  \-(El.k  "   El+l.k>      '    ^\E\i<n/2.k] 


with 


A|*Li/2.*  =    |A|(0i+lfJt  "  Qitk) 

Wi+i/2,k  =  MO)  ,0  =  QiQi+Lk  -  Qi.k) 
Flux  difference  splitting  can  be  accomplished  by  enforcing  interaction  of  neighboring 
states  using  Roe's  approximate  Reimann  solver,  where  the  interaction  is  enforced  through 
neighboring  cell  interfaces.  The  local  time  linearization  of  the  difference  vectors  AE1  for 
left  and  right  states  Ql  and  Qr,  respectively,  using  Roe's  method  is  obtained  as  follows 

A£*  =  A(0)*AC>     ,     A*  =  RA+R-1  (51) 

where  AQ  =  Qr  -  Ql,  a  =  (3£/30)  are  the  flux  Jacobian  matrices,  A*  is  the  eigenvalue 
diagonal  matrix  and  R,R_1  are  the  left  and  right  eigenvector  matrices.  The  intermediate 
state  0  is  calculated  using  the  Roe  variables  as 
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p  «  /pIp^ 


-=    "l/p!  +  unfP~R 


yfPl  +  >ft~R 


-  .  g^gL:  ^  (52) 


\[P~l  +  y/P~R 


The  numerical  integration  scheme  is  performed  using  a  fully-upwind,  factored,  implicit, 
third  order  accurate  numerical  scheme.  The  upwinding  is  perfomed  with  flux-difference 
splitting  using  Roe's  scheme.  For  complex  flow  fields  first  order  accurate  schemes  are 
too  diffusive  and  they  do  not  produce  solutions  with  reasonable  grid  spacings.  Here  a 
higher  order  scheme  is  used.   The  numerical  scheme  is  given  below: 

[j  + V*&+.*  +  4ffc.*>]p*  (53) 

iz  +  V*c3.*  +  Afc;,*  -  Jto-^A.J]*  x  (§ft  -  Qlk)   = 

-Re~1h^(Slk^1/2  -  Slk_1/2)p] 

H^  =  At/A£,  and  A,  V  and  8  are  the  forward,  backward  and  the  central  difference 
operators,  respectively.   The  quantities  Sul/2$k,  PUMn,  Slktl/2  are  numerical  fluxes. 
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The  accuracy  of  the  numerical  solution  is  improved  by  subiterations  to  convergence 
within  each  time  step  iteration.  The  approximation  to  0"'1  at  each  subiteration  is  Qp. 
When  p  >  2,  during  a  given  level  of  subiteration  to  convergence,  0P  =  0ntl;  but  when  p 
=  1  and  no  subiterations  are  performed,  then  Qp  =  £>n,  and  01"1  -  0n'1 .  By  subiterating  to 
convergence,  linearization  and  factorization  errors  can  be  eliminated  during  the  iteration 
process,  because  the  left  hand  side  of  Eq.(53)  can  be  driven  to  zero.  The  in  viscid  fluxes  £ 
and  f  are  evaluated  using  the  Roe  upwinding.  The  numerical  fluxes  for  a  third  order 
accurate  scheme  are  given  by 


1  - 

fi,k  =  h^i+i. k'Qi.k)    ~  -g  Idf  1+3/2. k  +  2df i+i/2. kl 


+  —  [2dfui/2tk  +  dfi-i/2.*] 


(54) 


or  in  terms  in  expanded  form 


Ei+in.k  =  ±tBlOi§k.tt.1/2,k)    -  5(ft-liJt,5i+i/a.jt>]  (55) 


-   ^[A£r(0ifJt,Oi+lfjb,(Si+1/a#Jt)    +  2AE-(QU2tk,Qi+lik,li+3/2tk)] 


+  ^  [2AE*  (Q^lfk,Qltk^1+1/2tk)    +  AeM^^Cm,,,^.^,,)] 


Using  the  notation  of  Eq.(54),  h  corresponds  to  the  first  order  accurate  flux,  and  the  df 
and  df  terms  represent  the  correction  terms  added  to  the  first  order  accurate  fluxes. 
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Higher  order  accurate  shock-capturing  schemes  have  some  limitations:  they  may  select 
a  nonphysical  solution,  they  produce  spurious  oscillations  and  they  may  develop  nonlinear 
instability  in  nonsmooth  and  discontinuous  flow  regions.  More  appropriate  shock- 
capturing  schemes  suitable  for  the  computation  of  weak  solutions  are  the  TVD  schemes, 
described  earlier.  For  this  simulation,  the  Osher-Chakravarthy  TVD  scheme  [Ref.  20]  was 
implemented.  This  scheme  uses  flux  limiters  to  impose  constaints  on  the  gradient  of  the 
fluxes.   The  flux-limited  values  are  computed  as  follows: 

dfi+3/2,k  =  minmod[df 2+3/2,  k'$dfi+1/2,k] 

dfl+uz.k  =  minmod[dfl+1/Zik,$dfl+2/2ik]  (56> 

dfi+\/2.k  =  minmod[dfl+1/2ik,$dfl-1/2,k\ 

df  1-1/2, k  -  minmod[dfi-1/2,k,$dfi+1/2,k] 
where  the  minmod  operator  is  defined  by 

minmod[x,y]  =  sign(x)  x  max[0,min{ fx /,ysign(x)}J  (57) 

and  the  comparison  parameter  (3  is  in  the  range  1  «s  (3  s  4.  Substitution  of  the  flux 
limited  values  given  by  Eq.(56)  into  Eq.(54)  makes  the  entire  scheme  TVD. 

The  flux  differences  AE  are  evaluated  by  the  Roe  scheme  described  earlier.   The 
viscous  fluxes  sitk.l/3  are  computed  using  central  differences  as  shown  below 
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Si.k+l/2         S[Ql,k+l/2>  (Q0  i.k+i/2'  ^1, k+l/2* 

Oitk+1/2    =    \(Qi.k   ~   Q**l.Jt>  <58> 

The  eigenvalues  of  the  flux  Jacobian  matrices  A  ,B  are  evaluated  as  suggested  in  Ref. 
24  to  prevent  expansion  shocks.    The  terms  A  ,C   of  the  left  hand  side  operators  of 
Eq.(53)  have  been  linearized  using  the  Steger- Warming  flux  vector  splitting.  A  validation 
of  the  viscous  flow  algorithm  for  a  transonic  airfoil  is  shown  in  Fig.  D2. 
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Figure  Dl .   Comparison  of  Turbulence  Model  Predictions  [Ref.  23] 
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Figure  D2 .   Viscous  Flow  Validation  [Ref.  23] 
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APPENDIX  E 
GRAPE 

The  following  is  a  summary  of  the  grid-generating  algorithm  GRAPE.  Reference 
25  provides  a  detailed  explanation  of  the  algorithm  theory  and  Ref.  26  is  a  detailed  user's 
manual.  GRAPE  is  an  algorithm  that  transforms  x  and  y  coordinates  in  real  space  to  | 
and  r|  coordinates  in  computational  space.  GRAPE  stands  for  GRids  about  Airfoils  using 
Poisson's  Equation.  r\  =  0  forms  the  inner  boundary  of  the  grid,  and  r|  =  r^  forms  the 
outer  boundary.  For  the  C-grids,  |  =  0  is  located  at  the  rear  boundary,  and  ^  moves 
forward  clockwise  around  the  airfoil  until  it  reaches  h,max  at  the  rear  boundary  again.  The 
mapping  functions  are  required  to  satisfy  the  following: 

£x*  +  %y  =  P  (la) 

1x*   +1yy  =  Q  (lb) 

The  following  relations  are  useful  in  transforming  the  equations  from  real  space  to 
computational  space: 

Sx  =  y„/J  (2a) 

i,  =  -VJ  (2b) 

nx  =  -yg/J  (2c) 

T]y  =  X|/J  (2d) 

where 

J  =  x^  -  x^y^  (2e) 
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Applying  Eqns.  (2)  to  (1)  yields  the  following  transformed  Poisson  Equations: 

ax^  -  2|ix4n  +  xx^  =  -J2(Px<:  +  Qx,,)  (3a) 

oyK  -  2jiysn  +  xym  =  -J2(Py^  +QyT1)  (3b) 

where 

<j  =  V  +  yn2  (3c) 

M-  =  x^  +  y^  (3d) 

T  =  X^2  +  y^2  (3e) 

Solving  Eqn.  (3)  for  a  particular  choice  of  inhomogeneous  terms  P  and  Q  (forcing 
functions)  and  a  particular  set  of  boundary  conditions  causes  the  grid  to  be  generated. 

There  is  latitude  that  exists  with  the  choices  of  P  and  Q.  Different  choices  of  P  and 
Q  result  in  different  grids.  In  an  effort  to  choose  P  and  Q  to  allow  reasonable  effort  in 
the  computation,  P  and  Q  are  redefined  in  terms  of  four  new  variables.  In  addition  to 
geometric  constraints  and  the  original  Poisson  equations,  six  equations  with  six  unknowns 
are  left.   These  equations  can  be  solved  in  a  straightforward  iterative  scheme. 

The  following  geometric  constraints  are  imposed  on  the  grid: 

1.  Spacing  along  £  =  constant  lines  between  the  body  at  r|  =  0  and  next  grid  node 
tj  =  1  is  specified  by  the  user 

2.  The  angle  of  intersection  between  the  body  and  £  =  constant  lines  is  specified  by 
the  user 

3.  Same  geometric  constraint  as  1)  but  applied  to  the  outer  boundary 

4.  Same  geometric  constraint  as  2)  but  applied  to  the  outer  boundary 

These  geometric  constraints  are  embedded  in  the  P  and  Q  terms. 
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The  program  uses  a  method  called  Coarse-Fine  Sequencing  (CFS).  This  method 
greatly  accelerates  the  numerical  convergence  of  the  program.  The  coarse  solution  solves 
the  grid  at  every  third  point  in  the  a  direction  and  every  third  point  in  the  6  direction. 
The  fine  solution  uses  a  cubic  spline  fit  to  find  the  rest  of  the  points. 

Reference  27  details  changes  made  to  the  original  GRAPE  code  to  allow  the 
generation  of  periodic  cascade  C-grids  for  turbomachinery  flow.  In  addition,  the  outer 
and  inner  boundary  routines  were  modified.  All  the  features  of  the  original  GRAPE  code 
were  retained. 
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APPENDIX  F 
GRAPE  INPUT  FILES 


$GPID1 

JMAX=19  9,KMAX=6  4,NTETYP=3,HAIRF=5,  JAIRF=3  16  ,  T1IBDST=3  , 

MDS=2,DSI=0.O01,JTEBOT=13,JTETOP=18  7,XLE=0,XTE=l,HOBSHP=4,NOBDST= 

=  1, 

XLEFT=-0.3, 

XRIGHT=1 . 5 , YBOTOM=-0 . 2  7  6 , YTOP=0 .276,  RCOR1(  =  0 .20, ALAMF=1 . 0 , ALAMR=0 

0. 

HORDA=4, 1, 

HAXITA=10, 10, JPRT=- 1,HOUT=0 

$END 

SGRID2 

XTFRAC=1 .  0  ,  ROTAIIG=0  .  0  ,  WAKEP=1 .  0 

SEtID 

SGRID3 

AIRFX  -  6.00000,5.99948,5.99891,5.99815,5.99667, 

5.99481,5.99264,5.99106,5.98939,5.98766, 

5.98500, 

5.982  00,5.97870,5.97507,5.97108,5.96668, 

5. 96185, 5. 95654, 5. 95069, 5. 94426, 5. 93719, 

5.92941,5.92085,5.91143,5.90108,5.88968, 

5. 87715, 5. 86337, 5. 84820, 5. 83152, 5. 81318, 

5. 79299, 5. 77 07 9, 5. 74637, 5. 71951, 5. 68996, 

5.65745,5.62170,5.58237,5.53911,5.49152, 

5. 43917, 5. 38159, 5. 3 18 25, 5. 25556, 5. 19286, 

5. 13016, 5. 067 47, 5. 00477, 4. 94208, 4. 87938, 

4. 81668, 4. 75399, 4. 69129, 4. 62860, 4. 56 59  0, 

4. 503 20, 4. 44051, 4. 37781, 4. 31512, 4. 25242, 

4.  18973, 4. 12703, 4. 06433, 4. 00164, 3. 93894, 

3.87625,3.81355,3.7508  5,3.68816,3.62546, 

3. 56277, 3. 50007, 3. 43738, 3. 37468,  3.  31198, 

3. 24929, 3. 18659, 3. 12390, 3. 06120, 2. 99850, 

2.93  581,2.87311,2.8104  2,2.7  4  77  2,2.68503, 

2.62233,2.55963,2.49694,2.43424,2.37155, 

2. 30885, 2. 24616, 2. 18346, 2. 12076, 2. 058 07, 

1.99537, 1.93268, 1.86998, 1.80728, 1.74459, 

1.68189, 1.61920, 1.55650, 1.49381, 1.43111, 

1. 36841, 1. 30572, 1. 24302, 1. 18033, 1. 11763, 

1.05  493,0.99  224,0.92954,0.86685,0.8  0415, 

0.74146, 

0.6787  6,0.61541,0.557  8  3,0.5054  8,0.4  57  89, 

0.41460,0.3753  0,0.33955,0.30704,0.27749, 

0.2506 3, 0.22621, 0.20401,0. 18383,0. 16548, 

0.14880,0. 13363, 0.1 1985, 0.107 32, 0.09 593, 

0.08557,0.07615,0.067  59,0.05981,0.0527  4, 

0.04  631,0.04  04  6,0.03515,0.03032,0.02592, 

0.0219  3,0.0183  0,0.01500, 

0.013  69,0.01112,0.009  87,0.00750,0.00536, 

0.00  3  51,0.002  01,0.00091,0.0002  3, 

0.00000,0.00023,0.00091,0.00201,0.00351, 

0.00536,0.00750,0.00987,0.01112,0.01369, 

0.01699,0.02062,0.024  61,0.02901,0.03  38  4  , 

0.03915,0.04500,0.05143,0.05850,0.06628, 

0.07483,0.08  426,0.09462,0.10601,0.11854, 

0.13232, 0.14749,0. 16417,0. 182 52, 0.20270, 

0.22490,0.24932,0.27887,0.31138,0.34713, 

0.38646,0.42972,0.47731,0.52966,0.58724, 

0.65058,0.71308,0.77608,0.83908,0.90208, 

0. 96508, 1.02508, 1.09 108, 1. 15108, 1.21408, 

1.27708, 1.34 008, 1.40308,1. 46608, 1. 52908, 

1.592 08, 1.652 08, 1.71508, 1.77808, 1.84108, 

1.904 08, 1.96708, 2. 03008, 2. 093 08, 2. 15 608, 
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2.21908,2.28208,2.34208,2.40508,2.468  08, 

2.53108,2.594  08,2.657  08,2.7  2  008,2.78  3  08, 

2.85000,2.92118,2.99378,3.05371,3.11665, 

3. 17956, 3. 24256, 3. 30553, 3. 36851, 3. 43149, 

3.49449,3.55748,3.61749,3.68049,3.74349, 

3. 80648, 3. 86948, 3. 93247, 3. 99546, 4. 05844, 

4. 12141, 4. 18137, 4. 24431, 4. 30725, 4. 371 08, 

4. 43309, 4. 49598, 4. 55886, 4. 62171, 4. 68455, 

4.74736,4.81015,4.87292,4.93566,4.99837, 

5.06105,5.12371,5.18633, 

5.24892,5.31849,5.38183, 

5.43945,5.49176,5.53935,5.58262,5.62195, 

5.65770,5.68899,5.71975,5.746  62,5.77104, 

5.79324,5.81350,5.83177,5.84845,5.86361, 

5.87740,5.88993,5.90132,5.91168,5.92109, 

5.92965,5.93743,5.94451,5.9509  4,5.95  67  8, 

5.9  6210,5.96741,5.97224,5.97664,5.98063, 

5.98426, 

5.98756,5.98930,5.99  097,5.99256,5.99474, 

5.99603,5.99810,5.99888,5.99945,6.00000, 

AIRFY  -  0.01500,0.01109,0.00940,0.00779,0.00557, 

0.00  3  65,0.002  09,0.0012  8,0.00066,0.0002  4, 

0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.000  00, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000, 0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.0000.0,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000,0.00000,0.00000, 

0.00000,0.00000,0.00000, 

0.00006,0.00051,0.00091,0.00201,0.003  51, 

0.0053  6,0.00750,0.00987,0.012  4  0,0.01500, 

0.01761,0.02013,0.022  50,0.024  64,0.0264  0,0.02799, 

0.02910,0.02940,0.02994, 

0.03015,0.03034,0.03061,0.03  089,0.03118, 

0.03150,0.03186,0.03  22  5,0.03272,0.03  316, 

0.03364,0.03426,0.03493,0.03559,0.03  639, 

0.03720,0.03813,0.03915,0.04025,0.04150, 
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$END 


0.04282,0.04436,0.04612,0.04814,0.05034, 

0.05272,0.05539,0.05830,0.06153,0.06502, 

0.06887,0.07273,0.07658,0.08043,0.08428, 

0.08814,0.09181,0.09584,0.09951,0.10337, 

0.1 07 22,0. 11107,0. 11493,0. 11878, 0.12263, 

0.12649,0. 13016, 0.13401, 0.13786,0. 14172, 

0.  14557,0. 14942,0. 15328,0. 15713,0. 1609 8, 

0.  16484, 0.16869,0. 17236, 0.17621,0. 18007, 

0. 18392, 0.18777, 0.19163, 0.19548,0. 19933, 

0.20342,0.20763,0.21145,0.21434,0.21709, 

0.21955,0.22172,0.22359,0.22517,0.22645, 

0.22745,0.22814,0.22854,0.22866,0.22850, 

0.22804,0.22728,0.22623,0.22489,0.22326, 

0.22133,0.21923,0.21673,0.21394,0.21085, 

0.2 07 47, 0.20380,0. 19984,0. 19558,0. 19103, 

0. 18619, 0.18106,0. 17563,0.16991,0.16390, 

0. 15760,0. 15101, 0.14412, 

0. 13695,0. 12859, 0.12073, 

0.11328,0. 10630, 0.09 977, 0.09367, 0.08802, 

0.08278,0.07811,0.07345,0.06931,0.06551, 

0.06201,0.05879,0.05585,0.05316,0.05064, 

0.04841,0.04634,0.04  444,0.04  271,0.04113, 

0.039  69,0.03781,0.03717,0.03608,0.03508, 

0.03417,0.03326,0.03242,0.03167,0.03098, 

0.03035, 

0.02  978,0.02937,0.02876,0.02796,0.02641, 

0.02451,0.02  230,0.02069,0.01901,0.01500, 
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SGRID1 

JMAX=1 

JAIRF= 

DSI=.0 

XLEFT 

$END 

SGRID2 

HOBCAS 

DSLE  = 

YSCL=1 

WAKEP= 

$EtID 

$GRID3 

AIRFX= 


69,  KMAX=40,  HTETYP=3 ,  HAIRF=5,  HIBDST 
316,  JTEBOT=  25,  JTETOP=145,  NORDA=0,3 
025,  XLE=0.0,  XTE=0. 61786, 
-0.3089,  XRIGHT=1.235,  RCORN=  0.0333, 


=0,  NLE=20,  NTE=  10,  DSRA=  0.5, 
0.0005,  DSTE=  0.0010,  PITCH=  0.50, 
.0  ,  XTFRAC=0.8  ,  ROTAHG=  -51.84, 
0.8, 


=7,  H0BSHP=7, 

MAXITA=   0,50,  N0UT=4, 


6. 

00000, 

5. 

99948, 

5. 

99891, 

5 

.99815 

, 

5.99667 

5. 

99481, 

5. 

99264, 

5. 

99106, 

5 

.98939 

i 

5.98766 

5. 

98500, 

5. 

98200, 

5. 

97870, 

5. 

97507, 

5. 

97108, 

5. 

96668, 

5. 

96185, 

5. 

95654, 

5. 

95069, 

5. 

94426, 

5. 

93719, 

5. 

92941, 

5. 

92085, 

5. 

91143, 

5. 

90108, 

5. 

88968, 

5. 

87715, 

5. 

86337, 

5. 

84820, 

5. 

83152, 

5. 

81318, 

5. 

79299, 

5. 

77079, 

5. 

74637, 

5. 

71951, 

5. 

68996, 

5. 

65745, 

5. 

62170, 

5. 

58237, 

5. 

53911, 

5. 

49152, 

5. 

43917, 

5. 

38159, 

5. 

31825, 

5. 

25556, 

5. 

19286, 

5. 

13016, 

5. 

06747, 

5. 

00477, 

4  . 

94208, 

4  . 

87938, 

4  . 

81668, 

4. 

75399, 

4. 

69129, 

4  . 

62860, 

4  . 

56590, 

4. 

50320, 

4. 

44051, 

4. 

37781, 

4. 

31512, 

4. 

25242, 

4. 

18973, 

4  . 

12703, 

4. 

06433, 

4  . 

00164, 

3. 

93894, 

3. 

87625, 

3. 

81355, 

3. 

75085, 

3. 

68816, 

3  . 

62546, 

3 

56277, 

3 

50007, 

3 

43738, 

3  . 

37468, 

3 

31198, 

3 

24929, 

3 

18659, 

3 

12390, 

3. 

06120, 

2 

99850, 

2 

93581, 

2 

87311, 

2 

81042, 

2. 

74772, 

2 

68503, 

2 

62233, 

2 

55963, 

2 

49694, 

2. 

43424, 

2 

37155, 

2 

30885, 

2 

24616, 

2 

18346 
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